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® A technique for the detection and removal of local defects in digital continuous-tone images. 
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© The present invention is a method for automatically detecting and correcting a wide range of local digital 
image defects with minimal user intervention. The detection process employs brightness and color thresholds in 
conjunction with magnitude thresholds on residuals of nonlinear spatial filters to separate defects from scene 
content with minimal confusion. The detected defects are then cosmetically corrected by combinations of 
nonlinear smoothing and grey-scale erosion. Several options are outlined for the feature selection, detection, and 
cleaning operations depending on source type and computational constraints. 
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Cross-reference to Related Applications 

The present application is related to U. S. Patent Application Serial No. 07/934,089 entitled "Process for 
Detecting and Mapping Dirt on The Surface of a Photographic Element" filed 08/21/92 by Robert Gray et 
5 al., corres. to EP 0 585 759 A1 . 

Fieid of the Invention 

The present invention is related to the field of scanning images, primarily from photographic film, to 
w form electrical equivalents of the scanned images and more particularly to the detection and the removal of 
defects due, for example, to dirt contamination or physical damage of the film. 

A portion of the disclosure of this patent document contains material which is subject to copyright 
protection. The copyright owner has no objection to the facsimile reproduction by any one of the patent 
disclosure, as it appears in the Patent and Trademark Office patent files or records, but otherwise reserves 
is all copyright rights whatsoever. 

Background of the Invention 

Digital images created from electronic scanning of continuous-tone photographic film often reveal 

20 visually objectionable defects due to dirt contamination or surface damage of the film. These defects 
usually must be digitally corrected by manual digital retouching (e.g. via "dustbusting" and other cloning 
techniques), which require that the operator visually identify each local defect in the image. 

A patent of interest for its teaching in this art is U.S. Patent No. 4.189,235. entitled "Test Device For 
Dynamically Measuring The Degree of Dirt Accumulation On Bank-Notes" by Guter et al. That patent 

25 describes a method for inspecting opaque web materials for dirt accumulation. A reflected signal from the 
web material is sensed by three adjacent photosensors as the material moves past the sensors. When the 
signal sensed by the center sensor is significantly different than from that sensed by the side sensors the 
sensed area is classified as containing dirt. 

In conventional optical printing systems, if dirt or a scratch appears on the photographic image, the 

30 photographic negative or slide may be cleaned and the print remade or the print itself may be retouched. In 
the case where the images are written to an optical disc, it is not convenient nor desirable to rewrite the 
image, because many images are written onto a disc before the images are viewed, and once an image is 
written it cannot be removed. Therefore, it is desirable to monitor the dirt present in photofinishing 
environments and to assess the effectiveness of any selected film cleaning methods. 

35 In U.S. Patent No. 4,907,156, entitled "Method and System for Enhancement and Detection of Abnormal 
Anatomic Regions in a Digital Image" by K. Doi et al. there is disclosed a technique for taking the 
difference between an "enhanced" and a "degraded" version of the input image to identify and remove 
structured anatomic background. A series of single thresholds are then applied to the difference images 
using the shape/size behavior of image regions to perform a feature extraction that is related to the 

40 detection of lung nodules. The feature extraction is implemented by the use of thresholds against which 
difference image pixels are compared to identify which super threshold "blobs" are lung nodules and which 
are not. 

A significant problem in the detection of dust is to provide an adequate way of detecting "small" 
objects. There are several possible approaches to this problem. For example, one could color-classify all 

45 the color pixels in an image, then spatially connect the pixels of a common color-class into segments, and 
then examine the spatial dimensions of the segment. This approach, while possible, has the drawback that 
the color classification must be quite accurate, or else some marginally-colored pixels will be misclassified. 
The problem is that differences in the direct RGB (or luminance-chrominance) image values do not alone 
give much contrast information between, say, a white small dust artifact and a bright grey extended scene 

so region on which the artifact is superimposed. Thus after color-classification one may be left with very large 
segments that are reasonably neutral and bright, but which now mask visible dust artifacts. 

A better approach prior to classification is to preselect only those points which have large local contrast 
from their surround. This may be directly accomplished by (1) spatially smoothing the color image and (2) 
subtracting the smoothed color image from the original color image. The resulting difference (or "residuals") 

55 image will emphasize all color pixels which have a large magnitude compared to some weighted local 
average magnitude. The key decision here is then the choice of spatial smoother. The most obvious first 
candidate is a low-pass linear filter (e.g. a simple spatially weighted averaging of all grey levels (GLs) in a 
local window); the properties may be easily analyzed mathematically, and (perhaps more important) such 
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filters may be rapidly implemented in software and hardware (the residual image is basically a high-pass- 
filtered version of the original). Unfortunately, linear residuals images have the property of retaining all local 
edges in an image, including very long edges due to scene discontinuities. The result is a very "cluttered" 
image (i.e. the residuals image contains many bright pixels having no connection with dirt or scratches); a 

5 significant amount of additional time must then be spent in scene reasoning to clean up the residuals map 
before defect decisions can be made. In short, linear-filter residuals poorly discriminate between local 
defects and true scene details. 

The identification of anomalies such as dust or scratches leads to the desirability of removing the same 
without the intervention of an operator. The present invention is directed to a methodology for identifying 

w and cosmetically correcting such anomalies with minimal operator intervention. 

Summary of The Invention 

One of the preferred methods of the present invention is: A method for the detection and removal of 
15 defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of pixel 
values expected, of defects; 

20 c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range of 
expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of residuals 
25 values expected of defects: 

f) forming a map in which each pixel that is within the range of expected defects according to steps b) 
and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

The present method, is called an autodustbuster (ADB) algorithm, and it has been successfully tested 
30 on both consumer-type color photographic images and color motion picture frames with good results. 

From the above it can be seen that it is a primary object of the present invention to provide an 
automatic local defect detection and removal technique for digital images. 

The above and other objects of the present invention will become more apparent when taken in 
conjunction with the following description and drawings wherein like characters indicate like parts and which 
35 drawings form a part of the present invention. 

Brief Description of The Drawings 

Figure 1 illustrates a system on which the present invention may be implemented. 
40 Figure 2 illustrates the windowing of an image raster to determine the grey level of a central pixel GLc 
at a position X. 

Figure 3 illustrates in chart form a histogram of GL values for an associated window. 
Figure 4 illustrates the averaging of the samples "nearest" in grey level GL to the central pixel GLc for 
the histogram of Figure 3. 

45 Figure 5 illustrates the replacement of the GLc value with the grey level average GL 3V at the position X. 
Figure 6 illustrates the process flow for a multi-band original image. 
Figure 7 illustrates the process flow for a resolution pyramid of an image. 
Figure 8 illustrates the process flow as applied to an R, G, B image. 
Figure 9 illustrates a flow for ADB processing within a resolution pyramid. 

50 

Detailed Description of the Invention 

Referring to Figure 1, the present invention is implemented using a digitizing scanner 10, a control and 
logic device 24, such as a computer with a display device 26, and an optical disc writer 36. As shown, a 
55 strip of film 15, containing one or more frames 18 of developed film, is placed into a scan gate 12 where 
the strip of film 15 is scanned under direction of control and logic circuitry 30. As each frame 18 is scanned 
the resultant scanned image data, represented by block 22, is digitized and transmitted to a memory 28 for 
storage. The computer 24 process the stored data, in a manner to be described, to provide output image 
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data 32 which may be written to an optical disc 34 by the optical disc writer 36 to provide a report as to the 
characteristics of the anomalies. The scanning device 10 is capable of quantizing pixel values into multiple 
brightness levels in separate red, green, and blue channels. A minimum number of brightness levels would 
be approximately 64 for an adequate quality image with a typical number being 256. Although the image 

5 data 22 has been described as originating from the operation of scanning film, it is well understood that 
other well-known techniques may be used to provide the image data. 

The autodustbuster (ADB) algorithm and several processing options that allow customization for 
different source types and computational constraints will now be described. 

The required ADB operations may be divided into (1) feature selection, (2) defect detection, and (3) 

10 defect correction. The operations must ensure that the maximum number of defects is detected and 
optimally corrected while minimally affecting scene elements which are not defects. Thus, the more 
confidence which the user has in the accuracy of the detection process, the more extreme can be the 
removal process; conversely, as one's confidence in avoiding "false alarms" from scene elements 
diminishes, milder removal methods may be used to avoid damaging real scene elements. 

is Visible dirt artifacts in digital positives, created from scanned color negatives, possesses the following 
properties: they are "white" (low color saturation), "bright" (very high apparent exposure values), "small" 
(small spatial extent) and "contrast" (large differences from their immediate surround). (If dirty color 
positives were instead scanned, then the dirt artifacts would be "dark" rather than "bright"). In addition, 
many such fragments exhibit little brightness variation within individual defects. 

20 It is obvious that the above characteristics are rather loosely defined, and that many real scene 
characteristics may possess one or more of these attributes, Thus the goals of an ADB algorithm are to 
maximize the percentage of defects detected ("hits") while minimizing the number of real scene elements 
mistakenly considered to be defects ("false alarms"). These goals tend to be competitive, that is, setting 
thresholds on features to increase the percentage of "hits" tends to increase the number of "false alarms" 

25 and vice versa. A detection mechanism must be used in which defects appear as dissimilar as possible 
from true scene elements. 

In most imagery, spatially small objects which differ greatly in brightness or color from the immediate 
surround are relatively rare. The ADB algorithm exploits this property by the innovative use of a particular 
form of nonlinear spatial filter. The properties of this filter are described in the next section. 

30 Figure 2, illustrates in flow chart form the generic autodustbuster process flow. The original image of 
block 100 is a multi-band digital image represented by one or more matrices of pixel values or equivalently 
grey levels (GL). In preprocessing block 110 some effective and convenient remapping of the original image 
may be performed. An example of preprocessing is the digital resampling of the image to reduce its size 
for computational expediency. An alternate preprocessing is the selection or recombination of the original 

35 image bands to produce a feature image with fewer bands. The detection suitability of the resulting feature 
image, in any case, should be at least comparable to that of the original image. Oval 120 represents the 
feature image that is outputted from block 110. The feature image is subjected to a band brightness test in 
block 130. The test applies low and high threshold evaluations to every pixel of each band in the feature 
image. Pixel values that are above the high threshold or below the low threshold in each band are 

40 considered to have passed the test and are forwarded to an edge-preserving-smoothing (EPS) filter 140. 
This filter non-linearly smoothes the selected pixel in all bands in a manner that retains major scene edges 
as subsequently described. The output of the filter at this pixel position is written to a corresponding 
position in the smoothed image 142. For those pixel values that are below the high threshold value, but 
above the low threshold value no subsequent EPS filtering is applied, i.e. the test is considered failed. Each 

45 band of the feature image pixel is copied to the resultant output of the smoothed image 142 and does not 
undergo EPS filtering. At difference node 145 the pixel values of each band of the smoothed image are 
subtracted from the pixel values at the corresponding position in the feature image 120 resulting in the 
residuals image 148. 

All of the pixel values from the residuals image are subjected to the residuals test of block 150. The 
so residuals test consists of a minimum and a maximum threshold that is applied to each band of the residuals 
image 148. When bright defects are to be detected, for a pixel to pass all band values must exceed the 
minimum threshold and at least one band value must exceed the maximum threshold. When a dark defect 
is to be detected a pixel will pass if all band values are less than the maximum threshold and at least one 
band value is less than the minimum threshold. 
55 Pixels that pass are forwarded to a classmap 160 and marked with a distinctive value. The classmap 
160 is a single band image representation which is used to mark the location of detected defects. Pixels 
that fail are marked with a different value in the classmap 160. In a classmap conditioning block 170 the 
values within the classmap 160 may be further modified to accommodate known spatial characteristics of 
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the defects. For example, a one-pixel morphological dilation of the defect-marked (passed) pixels in the 
classmap is typically performed. The conditioned classmap 175 is directed to the cleaning block 180 long 
with the original image pixel values. Pixels in the original image that are not marked as passed in the 
classmap will be copied unchanged to the output cleaned image 190. Pixels which are marked as passed, 
i.e., defect, are corrected using neighboring non-defect pixel values. The corrected values are then written 
to the cleaned image 190. 

A more detailed discussion of the major blocks of Figure 2 will now be undertaken. In ADB processing 
the feature image 120 used to detect the defect locations may be different from the original image 100 
which will be corrected. This is generally advantageous when computational time is at a premium and EPS 
filtering 140 of the full (generally three-color) image is too burdensome. In this case it is often possible to 
apply a preprocessing module 110 to create a one- or two-band "feature" image 120 which instead will be 
EPS-filtered and classified for defects. The effectiveness of this strategy will depend on the success with 
which a reduced-dimensionality feature image will retain the size, brightness, contrast, and color defect 
discrimination of the original image 100. This will of course be application-specific. For example, in digitized 
imagery originating from superimposed color cells (as in animated motion pictures) defects may be of any 
color and saturation; thus little is lost in detection performance using a reduced-dimension space jwhich 
sacrifices color and saturation information. In this case good success has been achieved in using a "band- 
maximum" one-band feature image 120 band where every feature pixel value is the maximum value of the 
input red, green, and blue bands. By contrast, in Photo CD applications the principle source of defects is 
dirt on the negative during scanning, which appears as "white" specks in the digital positive; here saturation 
information should be retained as an important identifier, and two- or three-band EPS-filtered feature images 
are the norm. Figure 7 illustrates the retention of all three color original color bands for classification. 

Rescaling of the feature image grey-level range during preprocessing 110 is often advantageous. In 
current implementations this consists of linearly scaling to a 0-255 pixel value range even for 10- or 12-bit 
original data. This results in a slight speed improvement during EPS-filtering 140 when the histogram is 
searched at every window location (see Figures 4 and 5 along with associated discussions). In Photo CD 
applications (e.g. Figure 7 this rescaling may be determined by finding .5% and 99.5% cumulative 
histogram points of the color-corrected bands of the entire input image. The minimum .5%-point and 
maximum 99.5%-point from among all bands are mapped to 0 and 255 pixel values, respectively in the 
feature image, thus ensuring that no color shift occurs. If the input data is not color balanced it may be 
necessary to perform a color balancing step on the original data. 

EPS Filters and k-Nearest Neighbor (kNN) 

To detect the critical size of the desired defects a smoothing filter is desired which retains extended 
image edges but smoothes small local contrasty image regions (blobs). Such filters exist as various 
nonlinear noise-reduction methods under the category of "edge-preserving smoothers" (EPS). The most 
widely known EPS filter is the median filter, which is one member of the class of rank-order filters. In rank- 
order filters, the pixel magnitudes are sorted and an average of some contiguous number of the sorted 
40 magnitudes is computed. For example, in a median filter the "middle" element of the sorted magnitudes is 
retained (an average over one sample). The sample mean is, at the other extreme, also a rank-order filter 
(although a poor EPS); it provides the average of all sorted magnitude values. 

A poor property of standard rank-order filters for image operations is that they are insensitive to the 
relative spatial location of the pixels within the current window; e.g. the magnitude chosen as the "median" 
45 may exist at the center of the window or at an extreme edge of the window. A number of specialized 
nonlinear window operators have-been devised to account for local spatial structure. One example is the 
sigma filter (J.-S. Lee, CVGIP.24,255-269. 1983). In this filter a histogram is created of the grey levels in a 
sliding window; based upon an estimate of the image noise, a range of histogram values about the grey 
level of the pixel at the center of the window is computed. One characteristic of the sigma filter, however, is 
so that it is specifically designed for noise removal and its averaging range is solely dependent on the noise 
characteristics; the filter averages over a different number of pixels at each window location. For "blob" 
detection, however, what is desired is a filter whose parameters are sensitive to local shape characteristics. 

A relative of the sigma filter which gives this control is the k-nearest-neighbor (kNN) filter (Davis and 
Rosenfeld, IEEE Trans, SMC, 7, 107-109 ,1978). This filter is implemented in the EPS filter module 140 of 
55 Figure 2. This filter has the following operation per color band: 

Referring to Figure 3, slide a rectangular window of odd dimensions n x m over an image in raster 
order. At each location determine the grey level of the center pixel (GL C ) or equivalents pixel value, and 
also create a histogram of the grey levels in the local window as per Figure 4. Beginning at the histogram 
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bin containing the center grey level GL C , sum contents of the histogram bins on both sides of the starting 
bin in the following fashion: Test if the. contents in the bin for GU are greater than or equal to k pixels. If 
this is true, then write out as the averaged value GL av just the existing value Gl^. tf the contents for bin GLc 
are less than k, then add the contents of the bin immediately below the bin and above GL^ to the sum and, 
again test whether the sum is greater than or equal to k. If the sum equals k, then calculate the average 
grey level according to the formula 

GL av = SUM {GLi X hist(i)} 



SUM {hist (i) } 



where the index i is over the summed histogram bins. If the sum is greater than k, then the extreme bins in 
the sum are proportionally reduced in weight in the calculation for GL av so that the total bin weight in the 
sum approximately equals k (as shown in the attached source code). If the sum is less than k ( then the next 
bin below GL^ and above GU are summmed and tested against k. This process is repeated until the sum of 
histogram bins about GL C are greater than or equal to k. Continue until k pixels (k < n x m) have been 
summed. Find the average grey level of the summed pixels per Figure 5 and write this value to an output 
image at the location of the window center per Figure 6. 

The k-nearest neighbor filter has some interesting filtering properties which are not immediately obvious 
from the above description. First, any constant-GL image structure (blob) which has an area of k or more 
pixels in the local window will not be significantly smoothed; any structure of less than k pixels will be 
increasingly smoothed as k decreases. Second, at "ramp" edges between otherwise constant grey level 
regions the slope of the ramp will be increased (sharpened) by the filter. 

The first property is the one of significance for "blob" detection. For example, if the window size is 5 x 
5 pixels, and if k is, say, 15, then extended straight edges having a depth (constant pixel value in both 
directions perpendicular to the edge) greater than three pixels will be effectively unchanged, and thus 
disappear in a residuals image. "Blobs" of area less than 15 pixels that fit entirely in the window will be 
smoothed and thus will appear strongly in the residuals image. This is the desired residuals property for a 
local defect detector, 

A major concern for a k-nearest neighbor filtering is computational speed. This is not a linear filter and it 
is not separable, nor does it have other usable symmetry properties. The major operations per window 
location are histogram creation and then accumulation with testing of the total accumulated samples. The 
current k-nearest neighbor implementation does significantly reduce the cost of computing the histogram by 
merely updating the window histogram in a horizontal strip of image rather than performing complete 
calculations per window. This results in a fractional saving of approximately (n-2)/n in histogram-computing 
operations (for a n pixel x m line window), i.e. the savings increase as window length increases. (This fast 
histogram updating is equivalent to that used in the "fast-median" filter-Huang and Tang, IEEE Trans, 
ASSP, 27, 13-18, 1979). 

Residuals Test 

As mentioned above, a residuals image 148 is created as shown in Figure 2. The anomaly classifier 
steps contained in modules 150 and 170 accept as input the residuals image 148 and create as output the 
final classmap 175 which marks the pixels to be cleaned. 

The residuals test in module 150 creates an initial classmap by performing tests on each pixel of the 
residuals image 148. In the classmap conditioning module 170 a spatial segmentation (i.e. connected- 
component labeling) of the non-background pixels in the classmap may be optionally performed with 
additional tests applied to the same segment. 

In the simplest (single-state) version and for detecting bright defects, each color residual pixel is 
evaluated against two residual-GL thresholds, GL rmin and GL rmax . In this implementation, all color residual 
pixels whose minimum band residual GL exceeds GL rrnin andwhose maximum band residual GL exceeds 
GL rmax are uniquely marked for subsequent cleaning; pixels which do not satisfy both of these constraints 
are marked to be left unaltered in the output image. (For single-band feature images GL rmin equals GUmax)- 
The two-dimensional map of category markings per pixel is termed the "class map". These thresholds 
impose a minimum required contrast between a defect pixel and its immediate background. The (smaller) 
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GL rmin threshold is generally set slightly above the "zero point" in the residuals image (e.g. to 133 for 
rescaled features with residuals biased by a value of 128 as described above) and acknowledges that valid 
bright defect pixels , e.g. due to dirt on a scanned negative, should not have a negative contrast (residual 
values less than the zero point) in any band. The (larger) GL^ threshold states a requirement of 

5 significant local contrast in at least one band for the defect to be visible. (Recalling that these thresholds are 
on the residuals GLs, it is obvious that not all bands should be required to pass a large GL rma „ value; for 
example, a white dirt artifact against a bright red background would have large residual values in the green 
and blue bands, but only a small residual value in the red band.) 

For two-state anomaly classification, two sets of residuals thresholds are used, i.e. GL rmin1 , GLrmaxi, 

10 GL rmin2 , and GL fmax2 - In general, the minimum-residuals thresholds are set equally, and slightly above the 
zero point as mentioned above; the maximum-residuals threshold is set somewhat higher for labeling high- 
confidence pixels than is the maximum-residuals threshold for the low-confidence pixels. This labeling 
results in a class map where the pixels possess one of three labels (normal, low-confidence dirt, and high- 
confidence dirt). In two-state classification, this class map undergoes image segmentation in module 170, 

is i.e. all residuals pixels spatially contiguous in the class map which have the same class are, assigned a 
unique code value. (These code values can be thought of as uniquely labeling all the "blobs" of suspected 
defects in the image.) A fairly efficient algorithm for segmenting the class map is employed which involves 
finding the lowest prototype code value in a spatial "tree" of prototype segments of the same class, and 
then renaming all the segments in the tree with this lowest code value. 

20 Following segmentation, all segments which contain low-confidence defect pixels and which are 
spatially touching segments containing high-confidence pixels are upgraded to high-confidence (i.e. the 
class map values of the pixels are reassigned to the "high confidence of an defect" class). The purpose of 
this. step is to make use of the adjacent presence of very obvious defects to decide that a slightly less 
contrasty region is part of the same defect. The segments which remain at low confidence are either 

25 relabeled as "clean" (i.e. all pixels of the segment are marked for no cleaning in the classmap) or else 
name their pixels marked for subsequent cleaning state in module 180. The segments, which are marked as 
high confidence of being defects, have their pixels marked in the classmap for subsequent erosion cleaning 
in module 180. 

30 Cleaning 

Two general cleaning methods are presently available for use within module 180 in the ADB algorithms: 
EPS substitution and greyscale erosion. In simple EPS substitution, the original input multi-band grey levels 
of image 100 are simply replaced by the band EPS-filtered grey levels of smoothed image 142 wherever a 
" 35 pixel is marked as a defect in the classmap 175 (as shown in Figure 2). This method has the advantage of 
computational simplicity in that EPS values have already been computed during defect detection; however it 
requires that the EPS image be stored for use during cleaning. In addition, EPS substitution is relatively 
forgiving of "false alarms", i.e. pixels classified as dirt artifacts but which are really part of the scene. 
Conversely, EPS-substitution cleaning may not completely remove all of an extended defect, to the degree 

40 that the defect itself influenced the local EPS values. 

In the greyscale erosion cleaning option, an iterative process occurs over detected defect regions, 
whereby GLs of defect pixels are replaced by a weighted average of original GLs of "normal" (i.e. non- 
defect or corrected defect) pixels within a small distance of the defect pixel. As the iterations progress, the 
effect of the "normal" GLs propagates toward the center of the original extended defect. The iterations end 

45 when all marked defect pixels have been replaced. Two precautions are in order. First, the "normal "-labeled 
pixels that are in immediate contact with the defect pixels at the beginning of the iterations are also 
relabeled to class "defect" prior to the onset of cleaning (via a morphological dilation operation which uses 
a 3x3-pixel probe function). This is because such pixels, while failing the detection thresholds, often have 
some contribution from the defect in their GLs which will be visible when propagated over multiple pixels 

so during erosion. Cleaning their values as well as those of "true" defect pixels can significantly improve the 
quality of the cleaning. The second precaution is that the spatial window over which the erosion weights are 
taken should be kept reasonably small, e.g. 5x5 samples. The reason for this is to minimize the influence 
of distant pixels which may not be part of the same scene structure as the region obscured by the central 
defect pixel. 

55 Current practice for the erosion cleaner is to use a 5x5-pixel sliding window with a minimum of six 
clean-class contained pixels in order for a central defect-class pixel to be filtered; if the central pixel is 
defect-class, but fewer than six clean samples are contained, then the defect classification is maintained to 
the next iteration. If the central pixel is of defect class and six or more clean-class samples are contained in 
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the window, then the GLs of the central pixel are replaced by a weighted average of the band GLs of the 
contained clean pixels, where the weights are linearly proportional to the inverse of the Euclidean spatial 
distance of the clean pixel from the central pixel. The cleaning is iterative (multiple passes through the class 
map) but is not recursive, i.e. the cleaning effect at one window location does not affect the cleaning at 

5 subsequent locations during the same iteration, in order to avoid directional biasing. 

Greyscale erosion has the advantage over EPS-substitution in being able to completely eliminate all 
trace of the artifact from the image; its disadvantages are increased complexity, a variable processing rate 
(the number of iterations is not known a priori), and the possibility for severe scene damage in the case of 
"false alarms" (i.e., real scene elements may be removed from the image). Despite these concerns, the 

70 excellent cleaning performance of erosion is currently the preferred method in both Photo CD (Figure 7) 
and Cineon (Figure 8) applications. 

Alternative Implementations: 

75 The ADB algorithm has several "flavors". These include (a) operation at a single spatial resolution, (b) 
operation within a resolution pyramid, (c) single-defect-state detection, (d) multiple-defect-state detection, 
and cleaning via (e) EPS substitution or (f) grey scale erosion. In the single-state versions, ail pixels of the 
image are classified as either "non-defect" or "defect", i.e. only one defect state is treated. In multi-state 
versions, the "defect" classification is further divided into subclassifications depending on the confidence 

20 that a pixel is indeed a defect; thus a three-state classification would consist of "non-defect", "high- 
confidence defect", and "transition-confidence defect". In single-state versions the cleaner used is either 
EPS-substitution or grey-scale erosion, but not both in one implementation. The multi-state version may 
employ a combination of both EPS-substitution for cleaning the lower-confidence defect regions and grey- 
scale erosion for cleaning the higher-confidence defect regions. 

25 The flavor combining options (a), (c), and (f) is a simple variant to implement, though it is not 
necessarily the fastest or best-performing version. In this case no spatial resolution pyramid is formed, the 
residuals image is not segmented, and the cleaner is chosen to be gray scale erosion for maximum 
correction. Figure 7 illustrates such a process flow and is typical of a Photo CD application. In this case the 
preprocessing 110 consists of rescaling the pixel values of the original image for maximum viewing contrast 

30 by linearly remapping the lowest .5-% band gray level to a value of zero and the highest 99.5% level to 
255, as previously described in the description for Figure 2. The brightness test 130 for dirt on scanned 
negatives consists of minimum threshold pixel values which must be equaled or exceeded in each band of 
the feature image 120 in order for the corresponding pixel to pass the test as a possible defect. The 
residuals classification 150 is a single state decision per pixel as either a non-defect or defect class, based 

35 upon a requirement that the smallest band residual equal or exceed the zero-background residuals level 
and that the largest band residual equal or exceed some user-determined threshold level which is above the 
zero background level. The classmap conditioning 170 in this case consists of a morphological dilation of 
the defect pixels in the classmap using a 3 x 3 or 5 x 5 kernel. The cleaning 180 consists of greyscale 
erosion as previously described under the section called "cleaning." 

40 Figure 8 illustrates an alternate processing flow which retains the characteristics, (a), (c) and (f), is 
appropriate for removal of bright defects in Cineon applications. In this case it is desirable to minimize the 
processing time required by the large image format the processing is therefor directed at only detecting the 
largest defects and at minimizing the number of bands in the feature image. The preprocessing 110 in this 
instance consists of creating a reduced resolution image rescaling it to a zero to 255 greylevel range such 

45 that the color balance is corrected, and then extracting a single band feature image in which the output 
pixel value is the larger of the Red, Green and Blue input values. The residuals test 150 applies a single 
threshold which each residuals pixel must equal or exceed to be marked as defect. Classmap conditioning 
170 consists of morphological dilation as described above. The cleaning module 180 consists initially of 
performing greyscale erosion of the reduced resolution unrescaled input image as previously described. 

50 Next the cleaned reduced resolution multi-band image is bilinearly interpolated back to the original size and 
the conditioned classmap is replicated to the original image size. Finally those pixels in the original image 
which are marked "defect" in the replicated classmap are replaced by the pixel values in the interpolated 
cleaned reduced resolution image. 

It is possible to use ADB within a multi-resolution pyramid structure to provide some speed gains from 

55 use of multiple smaller EPS kernels. This implementation would also allow more effective cleaning from the 
use of EPS residuals thresholds which are tuned to the specific pyramid level. 

Figure 9 illustrates a flow for ADB processing within a resolution pyramid. In most respects, ADB 
cleaning within a spatial level is similar to that in a single-level ADB version, i.e. detection, classification, 

8 



EP 0 624 848 A2 



and cleaning occur, and choices of one- or two-state and of EPS-substitution and/or greyscale erosion are 
the same. One key difference lies in the way the spatial pyramid residuals images are processed. 

In a spatial pyramid, an image is represented generally as a tow-resolution image plus a number of 
residuals images of increasing resolution. To reconstruct the image to a particular resolution, one interpo- 

5 lales the low-resolution image by a prescribed method to the next resolution level and then adds the 
residuals values of that level to each pixel GL. This process continues until the desired resolution level has 
been reconstructed. The Photo CD file format is an example of a spatial pyramid structure. Note that these 
residua! images are simply the pixel band difference between the original image at that resolution and the 
interpolated low-resolution version; they should not be confused with the EPS residual images discussed 

w above. 

In a pyramid structure, all higher-resolution reconstructed pixels depend on the lower-resolution GLs at 
that spatial location. Thus when some subset of defect pixels is cleaned at a given resolution level, the 
residuals values at all higher resolutions at these locations must also be adjusted; failure to correct the 
residuals' values could result in "ghost" images of the original defects upon pyramid reconstruction. 

75 Although various residuals processing schemes are possible, the method currently used for ADB is the 
following relatively simple one: the base (lowest level) version of the original image is cleaned as in Figure 
7. The classmap 175 prior to base cleaning is then replicated in module 200 to the size of the next highest 
resolution level. The defect-class pixels in the replicated classmap are then dilated in module 210 by a 
spatial extent equal to the region of support of the pyramid smoothing function for that level, i.e., the half- 

20 width of the kernal used to create the next-lower-resolution level. The resulting dilated classmap 220 is then 
compared with the pyramid-residuals image 230; the pyramid residuals image is assigned a zero pixel 
value in every band wherever the corresponding pixel in the class map 220 is labeled as defect. The 
resulting modified pyramid-residuals image 250 is then added per-pixel in module 270 to a version of the 
base level cleaned image 190 which has been spatially interpolated to the present pyramid resolution level 

25 in module 260. 

This results in a reassembled (base + 1)-level image 100' in which the results of the base-level ADB 
learning have been incorporated. This image is then ADB-processed starting with the preprocessing step 
110', and the process of Figure 9 thus repeats until all levels have been reassembled and cleaned. Note 
that ADB cleaning may be omitted from higher-resolution levels at the expense of small defects not being 

30 corrected, but that all higher-level, pyramid-residuals must be zeroed as described in order to ensure that 
no defect "ghosts" will appear in the final image. 

While there has been shown what are considered to be the preferred embodiments of the invention, it 
will be manifest that many changes and modifications may be made therein without departing from the 
essential spirit of the invention. It is intended, therefore, in the annexed claims, to cover all such changes 

35 and modifications as may fall within the true scope of the invention. 
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Scanner device 
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Scan gate 
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Film strip 
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Film cleaner 




18 


Frame 




22 


Image data 
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Control and logic device (computer) 
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Display device 
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Memory 
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Logic circuitry 




32 


Output image data 
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Optical disc 
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Disc writer 
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Original image 
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Preprocessing 
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Feature image 
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130 


Band brightness test 




140 


Edge-preserving-smoothing (EPS) filter 




142 


Smoothed image 
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Residual image 
Residuals test 
Classmap 
Conditioning block 
Conditioned classmap 
Cleaning block 
Cleaned image 
Classmap replication block 
Classmap dilation block 
Replicated dilated classmap 
Higher-resolution spatial resolution-residuals image 
Spatial residuals zeroing block 

Modified higher-resolution spatial resolution-residuals image 
Cleaned-image interpolation block 
Pixel value summing node 
Rebuilt modified higher-resolution image 
Higher-resolution preprocessing 

The Invention may be summarized as follows: 

20 

1 . A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
25 pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 

of expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

30 e) testing each residual value of step d) to determine if the residual value is within a range of 

residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

35 2. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

40 c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 

of expected defects; 

d) forming a residuals image which consists of the difference between the feature image pixel values 
and the filtered values of step c); 

e) testing each pixel value of the residuals image of step d) to determine if the residual value is within 
45 a range of residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

3. A method for the detection and removal of defects in digital images comprising the steps of: 
so a) preprocessing a digital image represented by pixel values to form a feature image in which defect 

pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values; 

55 d) forming a residual value for each pixel as a function of the difference between the corresponding 

feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 
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f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

4. The method wherein the map formed in step Q is conditioned by morphological dilation of the pixels 
s marked as defect. 

5. The method wherein the map formed in step f) is conditioned by morphological erosion of the pixels 
not marked as defect. 

6. The method wherein the map of step f) is formed with each pixel marked according to the magnitude 
of its values according to steps b) and e). 

w 7. The method wherein the map formed in step f) is formed with each pixel marked according to the 
magnitude of its values according to steps b) and e)and is conditioned by morphological dilation of the 
pixels marked as defect. 

8. The method wherein the map formed in step f) is formed with each pixel marked according to the 
magnitude of its values according to steps b) and e) and is conditioned by morphological erosion of the 

is pixels not marked as defect. 

9. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter. 

10. A method for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising the steps of: 

a) preprocessing a lower spatial resolution representation of the digital image by forming a feature 
20 image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

25 d) forming a residual value for each pixel as a function of the difference between the corresponding 

feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
30 b) and e) is marked as a defect; 

g) correcting the lower resolution representation of the digital image using the map created in step f); 

and 

h) correcting the higher resolution representations of the digital image using the map created in step 
f). 

35 11. The method wherein one or more of the higher resolution representations previously corrected using 
the lower resolution defect map of step f) are again processed using steps a) through h) to detect and 
correct for defects too small to be detected at lower spatial resolution representations. 

12. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter. 

13. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter. 
40 14. An apparatus for the detection and removal of defects in digital images comprising: 

means for preprocessing a digital image represented by pixel values to form a feature image in 
which defect pixels have large local contrast from their neighboring non-defect pixels; 

first means for testing the value of each feature image pixel value to determine if the value is within 
a range of pixel values expected of defects; 
45 filter means for edge-preserving spatial filtering of each of the feature image pixel values; 

means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the 'filtered value from said filter means; 

second means for testing each residual value from said forming means to determine if the residual 
value is within a range of residuals values expected of defects; 
so means for forming a map in which each pixel that is within the range of expected defects according 

to said first and said second testing means is marked as a defect; and 

means for correcting the digital image as a function of the formed map. 
15. An apparatus for the detection and removal of defects in digital images represented by multiple 
spatial resolutions ordered from lower to higher, comprising: 
55 means for preprocessing a lower spatial resolution representation of the digital image by forming a 

feature image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

first testing means for testing the value of each feature image pixel value to determine if the value is 
within a range of pixel values expected of defects; 
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filter means for edge-preserving spatial filtering of each of the feature image pixel values that are 
within the range of expected defects; 

residual means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the filtered value of said filter means; 

second testing means for testing each residual value from said residual means to determine if the 
residual value is within a range of residuals values expected of defects; 

mapping means for forming a map in which each pixel that is within the range of expected defects 
from said first and said second testing means is marked as a defect; 

means for correcting the lower resolution representation of the digital image using the map created 
by said mapping means; and 

means for correcting the higher resolution representations of the digital image using the map created 
by said mapping means. 

16. A method for the detection and removal of local defects in digital images comprising the steps of: 

a) generating a feature image from an original digital image; 

b) creating an EPS-residuals image using the feature image; 

c) testing and mapping image pixels based upon the values of the EPS residuals image and the pixel 
values of the feature image; 

d) modifying the mapped image as a function of mapped values of neighboring pixels; 

e) cleaning of defect pixels by an EPS-substitution process or by a greyscale erosion process; and 

f) cleaning of defect pixels in different spatial resolution representations of the original image if 
different spatial resolution representations are present. 
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Busteps 042193: (4-21-93) Applies a selective-averaging 

(x>NN smoothing) filter to a 3-band image; also . ■ , 

calculates the difference (residual image) between the original 
& the smoothed images - 

Search window size must be symmetric, odd of dimensions 

3x3 thru 51x51. ^ rt ^ . . 

Filter currently treats data as in range 0-4999 *oniy . 
Input image array is stored as Int"2. 

Intended for operation on a Sun workstation (e.g. Sparcstation 10) 
running Sun OS. 



Calls: 



75 



20 



25 



30 



35 



40 



45 



10 



1 

2 
3 
4 
5 



Read Hdr isl 

Timaae-f ormat-specif ic calls) 

Ihister 

(image-f crmat-specif ic calls) 

Ihismnmx 
Winbust_042193 

(image-f ormat-specif ic calls) 

Hist Init 

Hist^Jpdate 

Iselavh 
Add Image 

~ {image-format-specific calls) ^ = _ 

program Busteps_042193 
integer"2 npix,nlin 

integer*2 ibuf f (4096, 3) , buff n (4 096, 51 , 3) 
integer nav, npp, diagf lg 

Integer npixm, nlm, np, nl, band, nbands, nbx 
integer pixf ormi, pixf ormo, pixf ormro, duraml, dumm2 
integer gmin (3) , gmax (3) 
integer hist (5000, 3) 
integer*4 britelo, britehi 

real*4 *f actl, f act2, bias 

real*4 deltime, Dtirae, timearray (2) 

character*50 inf il , outf il, resf il 

data npm,nlm/51,51/ !Max EPS window size! m 

'data pixformro/3/ ! residuals image written as unsigned byte. 

data factl, fact2/l. 0,-1.0/ !Used in Add_Image! 

data npixm, nbx/ 4 09 6, 3/ 

data diagf lg/0/ !For winbust3! 

common /par/dumml, dumni2,nav 



format ( 



50 



Interactive parameter input 

write (6,5) 

-Busteps 042193-:'/, ^ ^ , , 

Tnput image may be 1,2, or 3 bands; /, 
Input values should be in range 0-4999; /, 
Max. image width: 4 096 pixels per line; /, 
Max. window height: 51 pixels'/. 
Left- i right -edge virtual mirroring /) 
write(6,10) 

formate ' Enter in f ilename [ch50] : ,5) 

read(5, ' (a50) ' ) inf il !Must be ISL IPP imagefile! 

write (6, 30) 



55 
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60 



70 



«»r EPS filename [ch50] [:' » s > 
' out resid u al s filename[ch50]:',S> 



30 
35 



50 

formatC 
51 read<5,«> nav 



read(5,*)np ^ae tfa , 

If tnP "write <«?•>' Window size out of rang 
Call Exit CD 

Endif 

nl=np 

npp = np m , , % 

write <6,S0> f avera ging points in window 

format (' I K 



"nils:-. '"S- jk'^? i— — 



stop 
■Endif 



write<6,60, fi high -brightnes S thresholds [in*].'." 

^i:-)*s"t.is:brit.h! 

^ ite ^! 7 Enter out residuals bias constant:',?) 
-^'•^islCinfii.npix.nXin^ands^ixfor.i, 

SS^^ifJ SSc * — "o iar.e-a.ort. 

stop 



Endif 



initialize arrays * constants 

npixm=npix 

nbx-nbands f d g/B/g2 , 

deltin.e=Dtime (timearray) 

pixformo=pixformi , Adtl ed 8/V92! 

pixforroro=pixfonni 

Calc histogram min.max of input image 

write <6, 100) infil Q) 

formate input file is prop erly filter the u*9- 

Tb is min/»ax call is *• required 

5ffi C ^|-«Sgi^^SfcS5SS-SSl«i- completer 

write ( 6 f * ) 1 

«i"(6^150) GLnin GLmax') 

150 format (/ ° 
Do ib«l.nbands 



100 
Note: 
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write {6, 151) ib, grain <ib) , gmax (ib) 
151 format <5x, 14, 2i 9) 

Enddo 

:..tesr if image mins & maxs are within LDT range of 0-4999: 
Do ib=l,nbands 

If (gmax(ib) .gt- 4999) then 

write (6,*)' Warning— gmax .gt. 4999 STOP 

Call Exit(l) 
Elseif (grain (ib) .It. 0) then 

write (6,*)' Warning— gmin .It. 0 STOP' 

Call Exit(l) 

Endif 
Enddo 

Note: If gmax gt 4999, then either rescale input data, or redimension 
his arrays in selav3 & selavb, & reset threshold on above tester. 



75 



Create EPS-filtered image 



20 



write (6,*)' [Beginning filtering]' 
If <np .ge. 3 .and. np -le. npm) then 

Call Winbust 042193 (infil, outfil, ibuff, buffn, . 

1 npixm, nEx,nlm, np, nl, ginin, gmax, pixf ormo, britelo, briteni, 

2 diagflg) 

Else 

Write (6,*) ' Np out of range stop' 

Endif 



Create EPS residuals image 

2 5 call Add Image (inf il, outf il, resf il, f actl, f act2, bias, pixf ormro) 

Tresiduals file is unsigned byte w/ mean^bias" ! 

write (6, 200) outf il 

200 formate Ending; EPS file is ',a50) 

write(6 f 201) resfil . ' , cm 

201 formate EPS residuals file is ',a50) 



30 



deltime=Dtime (timearray) 
write (6, 600) timearray (1) , timearray (2) 
80C format (/' User time (sec):', 

x f8.2,'; system time (sec) : ' , f 8 . 2/) 

end 



35 



Ibismnmx:-(5-28-91) Calcs. min & max values in a multiband histogram. 

Subroutine Ihismnmx (hist , nbands, gmin, gmax) 
40 Integer hist (5000, 3) , nbands, gmin (3) , gmax (3) 

Do ±b=l, nbands 

Do Jc=l,5000 

If (hist()c,ib) .gt. 0) then 
gmin(ib)=k -1 
go to 50 

45 Endif 

Enddo 

50 Do k«5000,l,-l 

If (Mst(k,ib) .gt. 0) then 
gmax (ib)«k -1 
go to 60 

50 



Endif 

Enddo 

60 continue 
55 Enddo 
return 
end 
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■ Read Hdr ISL: (9-25-90) Reads key info from header of an ISL-frat 
\ ~ image header- Closes image f iJ -f_^po^_f °_ lu .f i««-==== 

5 subroutine Read_Hdr_ISL (filename, npix, nlin, nbd,pxform) 

Integer imgnum, ios,npx,nln, nbd,pxform, rwmode 
IntegerM opnimg, getdef 
Integer*2 npix, nlin 
•Character*50 filename 

i _____ — — - — — __— _ 

jo . include ' /local/include/iopacxagef ' 

rwmode =0 

. . imgnum=0 _ nTCTr . 

ios=Opnimg (imgnum, filename, rwmode, -FALSE. ) 
If (ios -NE. SYSNRM) then ^ _ 

write (6,*) ' Error during OPNIMG— Return' 

75 return 
Endif 

ios=Getdef (imgnum, npx, nln, nbd,pxf orm) 
If (ios .NE. SYSNRM) then 

write(6,*)' Error during Getdef — Return 

return 

Endif 

Call Clsimg (imgnum) 
npix=npx 
nlin=nln 

return 
end 

25 
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• Ihister: (5-28-91) Calcs. multiband histogram of an image using a line 

• buffer only. Intended for use w/ seiectave2 . f or . 



8 



Subroutine Ihister (inf il , npixm, nbx, ibuf f , hist ) 
Character*50 infil 

5 integer arrtyp?ostatus, gstatus, rstatus, imgin, npx, nln, nbands, pixform 

integer opnimg, getdef , rdline, line 
integer hist (5000, 3) 
integer*2 ibuf f (npixm, nbx) 

: include ' /local/include/iopackage.inc' 

/0 * arrtyp=IDINT2 

-....Initialize the histogram: 
Do ib=l,nbx 

Do k=l,5000 

hist (k, ib)=0 

Enddo 

Enddo 

. . en the image: 
imgin~Q 

ostatus=opnimg (imgin, infil, 0, . false. ) 
If (ostatus .ne. SYSNRM) then 

write (6, 15) infil -*w«^-#x 
15 formate Ihister: input image ',a50,' not open—abort ) 

stop 

Endif 

gstatus=*getdef (imgin, npx, nln, nbands, pixf orm) 
If (gstatus .ne. SYSNRM) then 

write £6,20) _ . 

20 formate Ihister: Gstatus error on input — abort ) 

stop 
Endif 

:. ..Read the image lines into the histogram: 
Do j=l,nln 

30 rstatus^rdline (imgin, line, -1, ibuf f, npixm, arrtyp) 

If (rstatus -ne. SYSNRM) then 

write (6, 50) j - aH __,. M 

50 formate Ihister: error reading line ,i4, --abort j 

stop 

Endif 

35 Do ib=l, nbands 

Do i=l,npx 

)c=ibuff (i, ib) +1 

hist (k,ib)=hist ()c,ib) +1 

Enddo 

Enddo 

Enddo 
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Call clsimg (imgin) 

return 

end 
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job: winbust 042l93.v: 
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Winbust 042193: ^ ^ . , . _ . ^„ , 0 

Oses a barrel-shifting line buffer in multiple bands. 
Implements a "fast-histogram" method. 



70 



15 



20 



Key Pa f^f^ rS ^* 50/ name of input image file; not overwritten, 
outfil: Ch*50, name of output filtered image rile. 
0U int*2| all-band line buffer array; must be dimensioned by 

calling program . 
Int*2, strip buffer array; must be dimensioned by calling 

program. . . , -» ■ 

Int.* A, Max. allowed width of image, in pixels per line. 

Int*4! Actual width of image, in pixels per line. 
Read from image header. - 
Int*4, Max. allowed no. of bands in input image Inpuu. 
IntM, Max. allowed height of filter window (in lines) 
& thus max. allowed height of strip buffer. 

Int*4, Actual' width, height of filter window in pixels, 
lines. Input. , 
qm±n,gmax:IntM, 3-eleraent vectors. Assumed mm £ max band GLs 
gnu.n,g!n* , ^ t z in , age . Assumed previously computed. Input. 

Note this limits current max. no. of bands to 3. 
pixformo:Int*4, Datatype of output filtered IPP image file. 
pi ■ E.g. 3=uns byte, 4=signed Int*2, 6=signedR*4. 

Input . 

: IntM. Returned via Common. Set to current , , ac 

-band value of gmin & gmax, respectively. These values MUST 
be set prior to calling Iselavh; Function Iselavh uses 
these values passed via Common.. . 
IntM, no. of samples to average in each filter window. 
Passed thru Common by caller. 



ibuff : 

buffn: 

npixm: 

npix: 

nbx: 
nlm: 

np, nl : 



gmn, gmx : 



nav; 
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30 



General 
(1) 



(2) 



comments: , . . , _ ^ 

99% of the complexity of this subroutine is due to 
the implementation of the strip buffer. The advantage of the 
strip buffer is that it allows -convolution" rtype window 
operations while keeping at any one time in memor y only 
npix*nl*nbands worth of image, instead of npix*nlin*nbands for 
input & output images. . ne«»<; 

Note that the EPS filter is *N0T* recursive, i.e. it only uses 
CU>, original GL values in its input weights; . thus output from 
previous window filter do not affect neighboring window outputs. 
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Subroutine winbust 042193 (infil, outfil, ibuff , buff n,npixm, nbx, 
nlm, np, nl , gmin , gmax, pixf ormo, bntelo , br itehi , diagf lg ) 

Qiaracter*50 infil, outfil 

integer imgin , imgout , ost atus , gstatus 

integer nbx,nbands, pixf orm, pixf ormo noixb 

integer pst atus , o status2 , rstatus , wstatus , npx , nln , npixm, npix , npixc 

integer nlin,nlinb 

integer arrtyp, lininc, jtype 

integer gmn, gmx, Iselavh, nav, gmin (3) , gmax I J) 

integer nlm, np, nl , nl2 , np2, nbinx, gl_cent ^ 

integer opnimg, getdef , putdef , rdime, wrlme 

integer diagf lg 

integer*4 britelo, britehi, count 

integer*2 ibuf f {npixm, nbx) , buffn (npixm, nlm, nbx) 
integer*2 his (5000) 
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integers els (4096) !Note implicit max. of npix=4096 

data nbinx/5000/ 
5 common /par/cmn, gmx,nav 

include Vlocal/include/iopackage .inc' 
: rSve above comment-out before compilation! 

niS£/2 DINT2 »Half-height of filter window 
10 !E.g, for ni-5, nl2=2 ! 

n P 2-np/2 ! Half -width 

count=0 

•~~~toen~Input~fi output images for sequential I/O: 

75 ' * ' * 

imgin=0 

ostatus-SYSNRM , 
ostatus=opnimg (imgm, mi ii» .xaxse. > 
If (ostatus .ne. SYSNKM) then 



20 70 call Exit (1) 



?^(J 7 SLbSS"042 l 93:In image ...50.' not open-abort') 

CaJ 
Endif 

PtttuI=geS (imgin, npx, nln, nbands , pixf orm, 
If (gstatus .ne. SYSNRM) then 

write (6, 75) a hnrt') 

25 - K formate Gstatus error on input aoort ) 

Call Exit (1 J 

Endif 

npix-npx ! Actual no. pixels/line in image 

If .<npi*fc r :£. ( ^ Npixb > Npixm-Abort' 

Call Exit(l) 

Endif 

nlinb-nlin +2*nl2 

! nbands =Actual no. of bands in image 
80 format t/^ "Winbust 042193": Npix Nlin Type Nbands Typeo N P NaV , 

britelo,britehi,diagflg 
format (21x,2i5, i6, i7,i6, i4, i5,2i7, iB) 
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imgout=0 

pstatus^utdef Umgou^npx, nln, nbands, pi xfonno) 

It (Pstatus t .ne^S ^NKM^he n^ ^ ^.^^ 

Call Exit(l) 

45 Endif 

OStatUs2-SYSMFM ) 
ostatus2=o P nimg ^i^^f^rite-o n & 

If (ostatus2 -ne. SYSNRM) then ^ on . a hnrt' 
write (6,*)' Output image not open-abort 

50 
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Endif 



Call Exit(l) 



12 



lininc=nlin/8 iVariable used in writing current line-number 
jtype-lininc ! Ditto (just screen i/o related) 



rstatus=SYSNRM 
vstatus=SYSNRM 



70 

110 



.Load 1st strip into barrel buffer: 

Do j-l,nl2 +1 !Do per line in 1st buffer strip: 

line=j -1 ! (Damned -zero-starting" C array convention) 

rstatus=rdline (imgin, line, -1, ibuf f , npixm, arrtyp) 
If (rstatus .ne. SYSNRM) then 

write (6, 110) j , . , , 

formate "Winbust_042193" : Read error at line ',i4, '-abort ) 

Call Exit(l) 
Endif 



75 



! Transfer GLs from line buffer to strip buffer: 

! (Someone should be able speed up all this moving around of data) 
!i ~horiz pixel index for original image; t . 

!i2=horiz pixel index for central (non-mirror) pixels in strip bi 
!i3=horiz pixel index for edge (mirror) pixels in strip buffer. 

j2=j +nl2 
j3=nl2 +2 -j 
Do ib=l,nbands 

20 Do i=l,np2 +1 !Left 

i2=i +np2 
i3=np2 +2 -i 

buf fn (12, j2, ib) =ibuf f (i, ib) 
buffn (i3, j2,ib)=ibuf f (i, ib) 
buf fn (i2, j3, ib) =ibuf f (i , ib) 
buffn Ci3, j3, ib) *ibuf f (i, ib) 

25 Enddo 

Do i=np2 +2,npix -np2 -1 ! Central 

i2=i +np2 

buffn (i2, j2, ib) =ibuf f (i, ib) 
buffn (i2, j3, ib) =ibuf f (i, ib) 

Enddo 

. Do i=npix -np2,npix ! Right 

30 i2=i ^np2 

i3=2*npix +np2 -i 
buffn (i2, j2, ib)=ibuff (i, ib) 
buffn (i3 f j2, ib)=ibuff (i, ib) 
buf f n (i2, j3, ib) =ibuf f (i, ib) 
buffn (i3, j3, ib) =ibuf f (i, ib) 

Enddo 

35 Enddo 

! Write a status display to standard output! 

If <(j -nl2 -1) .ge. jtype) then el>Hift 
write(6,600) j !Write current line number to stdio. 

jtype= jtype +lininc 

Endif 

40 Enddo !End top - j- loop! 

jend=nl !"jend" is barrel index to most recent line 
jc«jend -nl2 !"jc" is barrel index to -middle" line 
Above sections only serves to fill the strip buffer for the 

first time; it gets more complicated as the buffer contents move 
"down" in the image. 
45 
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the buffer line index is rotated. 

Do j-nl2 + 2,nlin + n!2 .1 .!Indexes the r last- line +1) of current 

J j 0 =j - n i2 -1 ! Image index of output line! 

•Mark which pixels in the current line will be filtered: 
Do i=np2 +l,npixb -np2 

°° i0 - 1 ifS| f ;a.,c :is , :5:Sia s,-a& ,„...«.««. 

1 cls(i)=0 

go to 210 

Endif 

Enddo , , . , . _ , 

cls(i)=l !raark for filtering! 

coSnt-count +1 'Increment to-be-filtered pixel counter* 
0 continue 
Enddo 

Do ib=l, ribands 

gmn= grain (ib) 
gmx=gmax ( ib ) 

!Load window in central part of^barr el_buf£ff I 

i N 4te°thtt r ^xe?s-«i t hin i/2-winoow width of vert image boundaries 
!bc1& "-"ra'is above -i2". i.e. strip buffer horiz index. 
, , + , !Center-of -window index, 1st window! 

io£!-np2 ! Output pixel index! 

Call Hist init (nlm,np,nl,np2,buffn,npixa,nbx, 
1 npixb,nbands, i, jend,ib,nis,nbinx> 

g l_cent=buffn(i /j c : ib) i^o^windo^center 

If <ClS, iLfnioHbl=!selavh(his,nbinx.gl_cent) '.Filter 

Else 

... ibuff <io,ib)=gl_cent 

Endif 

!Update histogram for later window locations in a line: 

Do i=n?2 +2,npixb -n P 2 ! Cent er-of -window index! 
io=i -np2 'Output pixel index! 

Call Hist update(nlm,np,nl,np2,buffn, • 
x npixE, nbx, npixb, nbands, a , }end, ib, his, nbmx) 

!If passes brightness thresholds "St, then filter: 
gl_cent-buf fn (I, 3 c. ib) j^of^indow^enter 

If <C l DU ff(io?Ib}^ '.Filter 
Else 
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ibuff (io, ib)=gl cent 
Endif 14 



Enddo 



Enddo !End band loop! 

line=jo -1 

wstatus=wrline (iragout, line, -1, ibuff ,npixm, arrtyp) 

If (wstatus .ne. SYSNRM) then 

write(6,410) jo u 
formate w Winbust_042l 93" : write error at line ',i4, -abort ) 
Call Exit(l) 

Endif 

! Rotate barrel-buffer line indices: 

jend=jend +1 

If (jend .gt. nl) jend=l 

jc=jc +1 

If (jc .gt. nl) jc=l 

If (jo -le. nlin -nl2 -1) then 

j2=jend !For bottom mirror! 

■Read new image line: 
line=j -1 

rstatus-rdline (imgin, line, -1, ibuf f , npixm, arrtyp) 

If (rstatus .ne. SYSNRM) then 

write(6 r 430) j t . A t 

formate "Winbust 042193" :Error reading lxne ',i4,, -abort 
Call Exit (1) 

Endif 

Do ib=l,nbands ! Overwrite oldest line in barrel-buffer! 
Do i*l,np2 +1 'Left 
i2=i +np2 
i3=np2 +2 -i 

buffn(i2, jend, ib)=ibuf f <i,'ib) 
buf f n (i3 , jend, ib] =ibuf f (i, ib) 

Enddo 

Do i=np2 +2,npix -np2 -1 ! Central 

i2=i +np2 

buffn(i2, jend, ib)=ibuf f (i, ib) 

Enddo 

Do i=npix -np2,npix ! Right 

i2=i +np2 
13=2*npix +np2 -i 
buf f n (12, jend, ib) =ibuff (i, ib) 
buf fn (13 , jend, ib) =ibuf f (i, ib) 

Enddo 
Enddo 

Else ! Virtual mirror of bottom of image: 

j2=jend +1 'Buffer index of *top w line in buffer! 

If (j2 .gt. nl) j2=l 

j2=j2 -1 ! (12-18-92) 

If (j2 .It. 1) j2-nl 1(12-18-92) 
Do ib-1 , nbands 

Do i«l,npixb 

buf fn <i, jend, ib) -buf fn (i, j2, ib) 

Enddo 
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800 



Enddo 
Endif 

•Write a status display to standard output! 
if (j -9©- jtype) then 

ES&oS!* ^ne '.14.' is completed', 

jtype= jtype +lininc 

Endif 

Enddo iEnd main -j" loop! 



. Call clsimg(imgin) 

75 call clsimg (imgout) 

return 
end 
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17 

Hist init: (8-12-91) Used in Winbust_042193. Initialises * calcs a histogram 
" of a window. 

Subroutine Hist_init (nlra, np, nl,np2,buf f n, npixm,nbx, 
I npix,nbands,i, jend, ib, his, nbinx) 

integer*2 buf fh <npixm,nlm,nbx) 
integer*2 bis (nbinx) 

integer nlm, np, nl , n P 2 , npixm, npix, nbx, nbands , i , jend, ib, nbinx 
integer gmn,gmx, jO 

common /par/gmn, gmx 

...initialize histogram to zero: 
Do k*gmn+lf grox+1 
his(k)=0 

Enddo 

...Load histogram of window: 
j0= jend 

Do jj*l#nl 

10-jO +1 

If (10 .gt. nl) j0*l 
Do ii«l,np ' , 

i0*i -np2 -1 +ii 

ia-buffnUO, j0,ib) +1 

his(ia)«his(ia) +1 

Enddo 

Enddo 

return 
end 
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Job: hist update. vf 

Date: Wed 5pr 21 14:30:52 1993 



40 



45 



50 



55 



30 



Mcrwvm- -co 



EP 0 624 848 A2 



inteqer*2 buf fn <npixm,nlm,nbx) 
integer** bis (nbinx) 

, m nD nl np2 ,npixm,nbx,nbands,i, jend, ib,nbinx 
integer npm r nlm # np,nJL,r?p^#"F^ 
integer jO, idrop, iadd,npix 
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jO=jend 

idrop=i -np2 -1 
iadd=i +np2 

his(ia)«bis(ia) +1 

Enddo 

return 
end 
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: Iselavh: (6-19-91) Mod'd from selavb to eliminate heavy central 21 
; weighting on initial-point in interest of speed. 

; Should be valid for all window sizes. 

: Note that this function runs more slowly as the image dynamic range 

[gmax -gmin] increases. 
Assumes histogram is precalculated. 

"ac" is the GL of the "center" pixel in the current window. 
References: 

(1) "Noise Cleaning by Iterated Local Averaging", by L. Davis 
& A. Rosenfeld, IEEE Transactions Systems, Man, Cybernetics, 
Vol. SMC-8, No. 9, pp. 705-710, Sept. 1978. 

(2) w Digital image smoothing and the sigma filter", 

by J. S. Lee, CVGIP, V. 24, pp255-2€9 (1983). Sigma filter, however, 
uses a *fixed* range, in terms of +/- a*sigma, where sigma=est. 
std. dev. of additive noise. 

This function computes a selective-averaging filter output for a 
current location of a rectangular window of 1 band of an image array, 
& returns this output as the standard function return. 

Key Parameters: 

his: Int*2 (nbinx) ; Histogram of GLs in current window; 

Assumes that GL value corresponding to a histo bin is 
offset by. 1, e.g. his(l) corresponds to GL=0. 
nbinx: Int*4; Max. no. of bins in histogram, 
ac: Int*4; spatially-central GL in current window. 

gmn,gmx:int*4; min & max GLs in entire image for present band. 

Passed from caller via Common, 
nav: Int*4; desired no. of pixels to average in the 
current window. 

Output •iselavh" will be a weighted sum of histogram bin occupancies; 
this weighted sum is num/den. 



Integer Function Iselavh (his, nbinx, ac) 
real* 4 val, vail, val2 

30 

integer *2 his (nbinx) 
integer* 2 cnt 

integer gmin, gmax, nav, a c, nbinx 
Integer ar, al,denr, denl 
INTEGER* 4 NUM,DEN 

35 common /par /gmin, gmax, nav 

common /par2/num, den, ar, al, cnt, vail, val2, val 
i — — 



! average the nav-nearest grey levels . . . 
den-his (ac +1) 

if (den .ge. nav) then • !lf have >den local vals »ac,then no-change! 
40 iselavh=ac 

return 

endif 

nura~ac*his(ac +1) 
do k-1,4999 

ar-ac +)t ICurrent GL to "right" of center GL 

„ al-ac -k !- - "left" of - 

s if (ar .le. gmax) then 

num=num +ar*his (ar +1) 
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endif 
if (al 



endif 
if (den 



70 



den=den +his (ar +1) 

,ge. gmin) then 
num^num +al*his(al +1) 
den =den +his (al +1) 



.gt. nav) then 
cnt=den -nav 
denr=0 
denl=0 

If (ar .le. gmax) 
If (al .ge. gmin) 



22 
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else if 



Endif 



denr=his(ar +1) 
denl=his(al +1) 
val2-Float <denr +denl) 
vall=Float <ar*denr) +Float (al*denl) 
val=vall/va!2 

iselavh=Nint ( (Float (num) -cnt*val) /Float (nav) ) 
return 

(den .eq. nav) then 

iselavh=Nint (float (num) /float (den) ) 

return 



enddo 



20 



write (6, *) ' 
write (6,*)' 
write (6, *) ' 
write (6,*) ' 
write (6, *) ' 
Call Exit(l) 



Unexpected termination in iselavh!!' 
Local Iselavh values at Exit:' 
Central GL value (ac):' f ac 



Gmin, 
Nav: ' 



Gmax : ' , gmin , gmax 



nav 



end 



25 



30 



35 



40 



45 



50 



55 



34 



oncrwin- .CO ncojQjiaAi i - 



EP 0 624 848 A2 



23 



75 



20 



ggg g 


r rrr 


3 gg 


rr r 


? g 


r 


? g 


r 


3 gg 


r 


ggg g 


r 


g 




3 g 




gggg 






d 




d 




d 


aaaa 


ddd d 


a 


d dd 


aa; 


d d 


i a 


d d 


i aa 


d dd 


aaaa a 


ddd d 



aaaa 

a 

aaaaa 
a a 
a aa 

aaaa a 



d 
d 
d 

ddd d 
d dd 
d d 
d d 
d dd 

ddd d 



y 
y 
y 
y 

yy 
yyy y 
y 

y y 
yyyy 



ii 

i 

i 

i 

i 
iii 



m m nun 
mm ra in 
m m 
m m 
m m 
m m 



aaaa 

a 

aaaaa 
a a 
a aa 

aaaa a 



ggg g 
g gg 
g g 
g g 
g gg 

ggg g 
g 

g g 
gggg 



25 



30 



CCCC 



cccc 



aaaa 

a 

aaaaa 
a a 
a aa 

aaaa a 



t 
t 

ttttt 
t 
t 
r 

t t 
tt 



n nnn 
nn 
n 
n 
n 
n 



i 
i 
i 
i 
iii 



P PPP 
PP P 
P P 
P P 
PP P 
P PPP 
P 
P 
P 



35 



Job: add image. vf 

Date: Wed~Apr 21 14:28:56 1993 



40 



45 



50 



55 



35 



EP 0 624 848 A2 



24 

Add_Image: (6-26-91) Adds 2 1-3 band images w/ additive 
constants. 

Assumes the internal image representation is Real*4 
Max. # pixels/line=4096. 

Subroutine Add Image (infill, infil2, outfil, factl, fact2, const, pixformo) 
Character*50 infill, infil2, outfil 
Real*4 rbuf f 1 (4096, 3) , rbuf f2 (4096, 3) 
Real* 4 factl , fact2, const 

integer imginl, imgin2, imgout, nimx, ostatus, gstatus 
integer nbands, nbandsc, pixf orm, pixf ormc, pixformo, band 
integer pstatus, ostatus2, rstatus, wstatus,npx, nln, npx2,nln2 
integer opnimg, getdef , putdef , rdline, wrline 
integer arrtyp, lininc, j type, imax 

data imax/4096/ 

include ' /local/include/iopacJcagef ' 
ar rtype IDREA4 

Open input color image: 

imginl=0 

■ ostatus^opnimg (imginl , infill, 0, .false. ) 
If (ostatus .ne. SYSNRM) then 
write (6, 10) infill 
10 formate In color image ',a50, ' not open-abort') 

Call Exit(l) 

Endif 

gstatus »getdef (imginl , npx, nln, nbands, pixf orm) 
If (gstatus .ne. SYSNRM) then 
write (6, 20) 

20 formate Gstatus error on input — abort') 

Call Exit(l) 

Endif 



I f (npx . gt . imax ) then 

write (6,*) ' Npx .ge. imax — abort' 

write (6,*)' Npx=',npx 

write (6, * ) ' Imax= ' , imax 

Call Exit(l) 
Endif 
npix=npx 
nlin=nln 
write (6, 45) 

45 format (/' "Add Image": Npix Nlin Type Nbands Typeo' ) 

write (7, 50) npix, nlin, pixf orm, nbands, pixformo 
50 format (13x, 2i5, i6, i8, i8) 

...Open 2nd input class map: 
imgin2»0 

ostatus'opnimg (imgin2, inf il2, 0, .false. ) 
If (ostatus .ne. SYSNRM) then 
write (6, 60) infi!2 
60 formate In 2nd image ',a50,' not open-abort' ) 

Call Exit(l) 

Endif 

gstatus*getdef (imgin2 , npx2, nln2, nbandsc, pixf ormc) 
If (gstatus .ne. SYSNRM) then 
write (6, 65) 

65 formate Gstatus error on 2nd input — abort') 

Call Exit(l) 
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!? ^ki(6?M'« n ?nP^ r in.age 2 diSlAs?oSi 5K»1-C.U ExitU)' 

Call Exit (1) 
Endif 

pstatus^utdef Umgout , npx, nln, nbands, pixf ormo) 

If (pstatus .ne. SYSNRM) then _^™,r -abort' 

11 P write (€,*)' Pstatus error for output-abort 

Call Exit(l) 

ostatus2=opnimg (imgout, outf il, 1, .false. ) 
os i -i« for write-only! 

If <ostatus2 .ne. SYSNRM) then aHrt ^r ' 

11 * write (6,*)' Output image not open-abort 
Call Exit (1) 

Endif 

write(6,*)' Factl. Fcat2, Const=' , factl, fact2, const 

.Read each line of input image: 
lininc=nlin/4 
jtype^lininc 
Do i«l,nlin 

line«j-l 

IRead multiband line from color arr typ) 
rstatus^rdline (imginl, line, -1 , rbuf f 1, imax, arrxypj 

i for all hands! 

If (rstatus .ne. SYSNRM) then 

?o™Itr 8 Error reading image line -abort') 
Call Exit CD 
Endif 

!Read line from 1-band classniap: arrtvp) 
rstatus=rdline (imgin2, line, -1 , rbuf f 2, imax, arrtypi 

i --i" for all bands! 
If (rstatus .ne. SYSNRM) then 

15 ?o™Itr 8 !rror reading image line -abort') 

Call Exit(l) 
Endif 

Do band=l , nbands 

D ° 1= rbSffl(i band) -factl*rbuffl(i, band) 

rbufflU ' Dana, +fac t2*rbuff2(i,band) +const 

Enddo 

...Clip data to output range: 
If (pixformo .eq. IDBYTE) then 

Do i=1 ' n ?^ ffl(ifb and).toaxl(rbum 

rbuff 1 (i,band) «Aminl (rbuf f 1 (i, band) ,255.) 

Else n if° (pixformo .eq. IDINT1) then 

Do i=l^pix (i , band) mKmX (rbu ff 1 <i,band) , -128 ) 
rbuff 1 (^bandJ-Aminl (rbuffl (i,band) ,127.) 

Else^if 0 (pixformo .eq. IDINT2) then 

Do i-^ B P^ ffl(1#bim<l) . Anint (rbuffl <^M> 6fi v 
rbuffl <i,band) -Amaxl (rbuffl (i,band) , -32768 . ) 



0 



37 



EP 0 624 848 A2 



rbuf f 1 (i, band) =Aminl (rbuf f 1 (i, band) , 32767 . ) 26 

Enddo 

Else if (pixformo .eq. IDREA4 ) then 
Else 

vrite(6,100) ^ . 

5 10 o formate Fixf ormo-' , i3, ' not implemented-abort ) 

Call Exit (1) 
Endif 

Enddo !End band loop! 

c Write output line of data: 

line=j-l 

10 wstatus=Wrline (imgout , line, -1, rbuf f 1 , imax, arrtyp) 

If (wstatus -ne. SYSNRM) then 

write(6,*)' Error in wstatus-abort' 
. . Call Exit (1) 
Endif 

If (j .ge. jtype) then 

write{6,200) j ^ 
'5 200 format (10x,' Line ',i4,' is completed') 

jtype= jtype +lininc 

Endif 

Enddo !End line loop! 

z Close images 

on Call cisimg(imginl) 

20 Call clsimg (irngin2) 

Call clsimg (imgout) 

return 

end 
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Appendix B— Docket 66,719 

^ .or- N Y- 14650-2201 
© 1993 Eastman Koda* Company, Rochester, 
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c Bust Class Creates ** two-state ** class map from 3-band residuals image, 
c "Optionally dilates classmap before output. 

c Calls: Irgb_load_isl, Buster_Class_Create, Dilate_Class, Ioutput_isl. 

5 c Image arrays: 

c res ; 3-band dirty residuals image (input) 

c els ;l-band classmap (output) 

c residuals image & class map are assumed Integer*2 values (i.e., not rea: 

c Note demonstration max image size of 1B70 x 1374. 

10 c Does not assume that residuals image is biased by a predetermined 

c constant; this info is implicit in the GL-res thresholds entered by the 

c Calls: 

c Irgb load isl 

c ~" (calls image-format-specific io routines) 

c Buster Class_Create 

75 c Dilate Class 

c IoutpuE isl . . , • » 

c "(calls image-format-specific io routines^. 

Program Bust_Class 

lnteger*2 res (3, 1870, 1374) 
Integer*2 els (1870, 1374) , cls2 (1870, 1374) 

integer*2 npix,nlin 

Integer npixm, nlinm f band, nbands , nbc, nbx, cnt , np # nl 
integer pixf ormi,pixf ormo 

integer cnttot 
integer ixO, iy0,nx r ny 

integer dgl, glmin, glmax, gltmin, gltmax, dilf lg, byte 

real* 4 deltime,Dtime, timearray (2) 

character*50 inf il, outf il, resf il, clsf il 
30 character* 1 iname (50) , iname2 (50) 

character*50 commSO 

equivalence (commSO, iname2 (1) ) 

data dgl/255/ JDefines "dirt" class! 

data byte/3/ 

data nbx, npixm, nl iron/ 3, 1870, 1374/ 
35 data pixformo/3/ ! Classmap of type "byte" 

I — — 

write (6,*)' [Max. image size: 1870 x 1374]' 
write (6,*)' ' 
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write (6, 20) , ^ , _ . _ n . , 

20 formate Enter in biased residuals image filename [chbOj. 

read (5, ' (a50) ' ) resfil 

write (6, 30) , _ 

30 formate Enter out classmap image name:',$) 

read<5,' (a50) ' ) clsfil 

write(6,50) l v , f , 

50 formate Enter min, max residuals GL thresholds:',?) 

45 read (5,*) gltmin, gltmax 

write (6, 60) 
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2 

60 * formate Do you wish prefilter dilation? [yes= w l") : ' , $) 

read (5,*) dilflg 

5 C write "(6,*) 'ciass threshold used is ',dgi 

deltime=Dtime (time array) 



70 



!Read biased color-residuals image: 

Call Irgb_ load isl (resfil, res, nbx, nbands, npixm, nlinm, 
1 npix,nTin,pixformi, 0) 

write (6, 150) resfil 
150 formate Input residuals file is ',a50) 



t 



75 'Compute class map: 

Call Buster_Class_Create (res, cis, nbx, nbands, 
1 npixm, nlinm, noix, nlin, cltmin, gltmax, dgl) 



i 
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! Optionally, dilate classmap by one pixel: 

If (dilflg -eq. 1) then 

write <6,*)' [Dilating class map]' 

Call Dilate_Class (els, cls2, 
1 npixm, nlinm, npix, nlin, dgl) 

Call loutput_isl (clsfil, cls2, npixm, nlinm, 
1 npix, nlin, byte, iname) 

write(6,200) clsfil 
200 formate Out dilated class file is ',a50) 

30 Else ! write undilated classmap! 

Call loutput_isl (clsfil, els, npixm, nlinm, 
1 npix, nlin, byte, iname) 

write (6, 300) clsfil 
300 formate Out undilated class file is ',a50) 

35 Endif 

deltime=Dtime (timearray ) 
write (6, 800) timearray ( 1) , timearray (2) 
800 format (/' Bust Class:: user time (sec):', 
1 f 8. 2, '7 system time (sec) : ' , f 8 . 2/) 
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end 
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Job: irgb load isl dugas.vf 
Date: Wed Xpr 2T 14743:48 1993 
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c Irgb_ load_isl: Loads up to a 3-band image into a Int*2 array. 

Subroutine Irgb load_isl (bname,b,nbx,nb, npixm, nlinm, npix,nlin, pixform, 

1 lunprt) 

Character*50bnarae, infil 

Integer*2 b (nbx, npixm, nlinm) ,npix, nlin, ibuf f (4096,3) 
integer imgin, gstatus , ostatus 
integer nbands, pixf orm, band, nbx, nb, npixm, nlinm 
integer rstatus , npx , nln 
integer opnimg , getdef , rdl ine 
integer arrtyp, lininc, jtype,pixdim 

Data pixdim/4096/ 
include ' /local/include/iopackagef ' 

arrtyp=IDINT2 

c. . .Open input image : 
15 imgin*=0 

infil=bname 

ostatus~opnimg (imgin, infil, 0, . false . ) 
If (ostatus .ne. SYSNRM) then 
write (6, 10) infil 
10 formate In image ',a50,' not open-abort' ) 

stop 

so Endif 

gstatus«getdef (imgin, npx, nln, nbands, pixform) 
If (gstatus .ne. SYSNRM) then 
write<6,20> 

20 formate Gstatus error on input — abort ) 

stop 

25 Endif 

If (nbands .gt. 3 -or. nbands .gt. nbx) then 

write (6,*)' Input no. bands too large — abort' 
stop 

If (npx .gt. npixm .or. npx .gt. pixdim) then 

write(6,*)' "IRGB_load_isl" : too many pixels — abort 
30 stop . ~ 

Endif 

npix=npx 
nlin=nln 
nb**nbands 
write<6,45) npix, nlin, nbands, pixf orm 
35 45 format(2x,' "Irgb load isl":: Npix:',i5,'; Nlin : ' , i5, ' ; Nb: ,i2, 

1 Type:'7i2) " 

If (lunprt .gt. 0) then 

write (lunprt, 45) npix, nlin, nbands, pixf orm 

Endif 

c Read each line of input image: 

40 linii*c«nlin/4 

jtype-lininc 
Do j~l,nlin 

c ... Read input lines & load 3-d array: 

line=j-l . jt % . 

rstatus=rdline (imgin, line, -1, ibuff ,pxxdim, arrtyp) 
i f or a ix bands! 

45 If (rstatus .ne. SYSNRM) then 

write (7, 80) j 

BO formate Error reading line -abort') 

stop 
Endif 
Do i«l,npix 
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Do band=l,nbands 

b (band, i, j) =ibuf f (i, band) 

Enddo 
Enddo 

If (j ,ge. jtype) then 
vrite(6,200) j 

200 format (10x,' Line ',±4,' has been read') 

jtype-jtype +lininc 

Endif 

Enddo !£nd line loop! 

c. . .Close image: 

Call clsimg (imgin) 

return 

end 
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Job: buster class create dugas.vf 
Date: Wed Apr 21 1T:46:32~1993 
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&2£p e 2e»lS patently only bistate (dirt or ™ n ££ enc * ut 
may be generalized to multiple-conf ldence states of presence 
of defects. 

Subroutine Buster Class_Create (res, els, nbx, nbands, npixm, nlinm, 
npix, nlinTgltmin, gltmax, dgl ) 

Integer *2 res (nbx, npixm, nlirun) , els (npixm, nlinm) 
Integer*2 npix r nlin 

Integer nbx r nbands, npixm, nlinm, gltmin, gltmax, dgl, glmin, glmax 
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foroatr^ust Class_Create: Gltmin Gltmax Dgl') 



write (6, 11) gXtmin, gltmax, dgl 
* c# 217, i5/) 

Do j«l,nlin 



format (19x, 



Enddo 

return 

end 



Do i«l,npxx t 

glmin«res(l,i,3) 
glmax«glmin 

Do ib«2, nbands , J ... . . 

glmin-MinO (glmin, res (ib, i, j) ) 
glmax»MaxO (glmax, res (ib, i, 3) ) 

?f dC (glmin .ge. gltmin .and. glmax .ge. gltmax) then 

cls(i, j)-dgl 
Else j ^ 

cls(i, j)=0 

Endif 

Enddo 
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Job: dilate class dugas.vf 
Date: Wed Apr 21 lT:47:31 1993 
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c Dilate Class: Dilates by 1 pixel *" 1 f«« d 3 S a " ^'window, 
c Byte class-map. Dilation kernel is a 3x3 square 
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ctar-class which is dilated; if current pixel is undilated then 
old classvalue is copied as new classvalue. 

A1S ° ^SSol^^^S^ window still covering image. 



TSroutine Di latere lis s (clsi, clso, npixm, nlinm, npix, nlin, ctar) 

Integer *2 clsi (npixm, nlinm) , clso (npixm, nlinm) , npix, nlin 

Integer npixm, nlinm, ctar, im, ip/ jm, jp 

i 

ro™Itt/'°Di!ate_Class: target class is -.15/) 

Do j«l,nlin 

jm« j-1 

jm«Max0(jra, 1) 
jP=j+l 

jp~Min0 (jp, nlin) 
Do i«l,npix 

im«i-l 

im«Max0 (im, 1) 
ip=i+l 

ip»MinO (ip, npix) 
Do jj«jm,jp 
Do ii=im, ip 

'If any neighbor^ctar, then 
! set center=ctar as well: 

If (clsi<ii,jj) -«q- ctar) then 
clso<i,1)*ctar 
go to 100 

Endif 

Enddo 
Enddo 

!No change, so copy to output: 
clso(i, j)«clsi(i, j) 
continue 



Enddo 



return 
end 



Enddo 
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Job: ioutput isl dugas.vf 
Date: Wed Apr~21 T4: 48:35 1993 
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11 

r SUBROUTINE: IOUl-rUT isl (10-9-90) 

C PROGRAMMER: ELIZABETH J. DONALDSON/R- Gray 

C DATE: 9/20/90 

C OFILE - INPUT - CHAR*50 - DESTINATION IMAGE 

r FOUT - INPUT - Integer*2 - DESTINATION ARRAY 

C NPIXM - INPUT - INTEGER* 4 r MAX PIXEL DIMENSIONS FOR FIN 

C NLXhS - INPU? - INTEGER* 4 - MAX LINE ?™ENSION FOR FIN 

C NPIX - OUTPUT - INTEGER* 4 - NUMBER PIXELS IN JNTILE 

C NLIM - OUTPUT - INTEGER* 4 - NUMBER LINES IN INFILE 

C PIXFORM - OUTPUT - INTEGER* 4 - INFILE DATATYPE 

C INAME50 - OUTPUT - CHARACTER UNUSED^ 

inameSO) 

include ' / local/include/iopaclcage . inc' 



Subroutine Ioutput_isl (of ile, f out, npixra, nlinm, npix, nlin, pixf orm, 



character of ile*50, inameSO (50) 

75 lnteger*2 f out (npixm, nlinm) 

integer*2 npix, nlin 

integer putdef , opnirng, wrband . - < 

integer*4 imgout , arrtyp, gstatus, pixf orm, npixm, nlinm 

integer*4 ostatus, pstatus, wstatus, npx, nln 

i _ 

write (6,5) npix, nlin, pixf orm , , 

20 5 formats loutput isl: image header data: np « nx f ±*, 

+ ' pixf orm * ' , i3T 

arrtyp ~ IDINT2 
imgout=0 

C DEFINE CHARACTERISTICS OF OUTPUT IMAGE 

25 npx=npix 

nln=nlin 

ostatus = putdef (imgout, npx, nln, 1, pixf orm) 
if ( gstatus .ne. SYSNRM) go to 1000 

C OPEN OUTPUT IMAGE FOR RANDOM ACCESS 

ostatus » opnimg (imgout , of ile, 1 , -true . ) 
30 if ( ostatus .ne- SYSNRM) then 

write (6, 10)ofile 
10 "formate Output image' , a50, ' not open. ) 

write (6,*) ' Exiting — ostatus= ostatus 
go to 1000 

endif 

35 if ( pixform.eq.IDBYTE) then 

do 100 j = l,nlin 

do 200 i - l,npix 

if (foutU, j).lt.0) fout(i,j) = 0 
if (foutU, j) .gt.255) fout(i,j) = 255 
200 continue 
100 continue 
40 elseif( pixf orm. eq.IDINTl) then 

do 110 j - l,nlin 

do 210 i - l,npix 

if (foutU, j) .It. -128) fout(i,j) ■ ~i 
if (fout <i, j) -gt.127) fout(i,j) - 127 
210 continue 
110 continue 
45 elseif( pixform.eq.IDINT2)then 

do 120 j - l,nlin 

d °il 2 ( fo«(iS-lt.-32768) /out(i j) --32768 
if ( fout <i,j) .gt. 32767) fout (i,j) =32767 

50 
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220 continue 
120 continue 

elseif ( pixform.eq.IDREA4) then 
GO TO 900 

else 

s write (6, 30)pixform 

30 format <' IOUTPUT ISL: Pixform ',i3, ' not implemented.') 

endif 

900 wstatus « wrband<imgout, 0, f out,npixm, nlinm, arrtyp) 
if ( rstatus .ne. SYSNRM) go to 1000 

10 1000 continue 

call clsirog (i ragout) 

return 

end 
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Appendix C--Docket 66,719 
© 1993 Eastman Kodak Company, Rochester, N.Y. 14650-2201 
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Job: buster erode am dugas.vf 
Date: Wed Apr 21 l5:lT:14 1993 
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: Buster Erode am: (11-11-92) Uses dynamic memory allocation on the Sun. 

: Calls Hedirt5x5 byte2 nopad; no 2-pixel -dead zone" 

; around image; aTl pixels in image are potentially cleaned. 

; Intended for use on a Sun workstation running SunOS (e.g. a 

; Sparcstation 10) . 

5x5 erosion window is currently hard-wired in. Minimum clean fraction 
of erosion window is nominally .24. 

Cleans RGB image by grayscale erosion 
using an external 2-state class map. 

Calls: 

Irgb load isl 

~~ (includes image-f ormat-specif ic calls) 
Binput Isl 

~" (includes image-format-specific calls) 
Dedirt5x5_byte2_nopad 
Irgb write_isl 

(includes image-format-specific calls) . 

Image arrays: 

iml : Original (dirty) 3-band image (input) 

els :l-band classmap (input) 

\rfi7 : Cleaned 3-band image (output) 

Warning: Once erosion iteration begins, iml & im2 contents 
are *partially*-cleaned image; original is not retained. 

Removes by interpolation all pixels whose class-map GL is .ge. 
a threshold value. Works on a multiband image (up to 3 bands) . 
Assumes iterative in/of recursive filtering. 

Uses adjustable variable-input window size, allows only filtering 
a rectangular piece of the image . * 
Image & class map are assumed Integer*2 values not re f*ii_ 



30 



35 



Program Buster_Erode_Am 

Pointer <piml,iml), (pim2,im2), (pels, els), (pcls2,cls2) 

Integer*2 iml (1, 1, 1) , im2 (1, 1, 1) 
Byte cls(l,l),cls2(l,l) 

Integer*2 iml (3, 1828, 1332) , im2 (3, 1828, 1332) 
Byte els (1828, 1332) , cls2 (1828, 1332) 

integer*2 npix,nlin 

Integer npixm, nlinm, band, nbands, nbc, nbx, cnt 
integer pixf ormi, pixformo,pixf ormc 

integer iml_bytes, im2_bytes, cls_bytes, cls2_byt'es, tot_bytes 
integer itmax, it, cnttot 
integer inflg 
integer dgl,byte 
Integer*4 iystart, iyend 

real*4 deltime, Dtime, timearray (2) , fmin 

45 character*50 infil, outf il, clsfil 

data itmax/40/ !Max. no. of erosion iterations! 

data dgl/255/ 'Defines "dirt* class! 

data diagflg/0/ 
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data byte/3/ 

data nbx,npixm, nlinm/3, 1828, 1332/ 



10 



write (6,*)' 77 , 

write (€,*)' Buster Erode_Am : Dynamically allocated memory 

write(6,10) , wcA1 # e , 

formate /' Enter in "dirty" file nametchSOj: ,5) 
read (5,' (a50)')infil 



Call Read hdr isl (inf il,npix,nlin,nbands,pixf ormi) 
write<6, 15) npix,nlin, pixformi, nbands 
15 formate npix-',i4,'; nlin«',i4,'; type=',12,'; bands= ,12) 



npixm=npix 
nlinm=nlin 
nbx=nbands 



irol byteso2*nbx*npixn*nlinm 
im2~"byt e s= iml_byt e s 
c 1 s~by t e s«np i xm * nl iron 

cls7 bytes=cls bytes , w a , 0 . _ ^ 

tot Eytes«=iml 5ytes +im2_bytes +cls_bytes +cls2 bytes 
write (6,*) '"Memory allocation requested (bytes): ',tot_bytes 

write (6, 20) 

20 formate Enter in class filename:' r $) 

read (5,' <a50) ' ) clsfil 

write (6, 30) wcft , 
30 format ( ' Enter out "cleaned" f ilename [ch50] : ,9) 

read<5 f ' <a50) ' ) outf il 

- write (6, 40) _^ , . ._ , c . 

40 formate Enter reqd cleaning fraction -fmin" [Real]: , $) 

read (5,*) fmin „ C1 . 

If (fmin -le. 0. .or. fmin .gt. 0.5) then 

write (6,*)' Fmin out of range — abort 

Call Exit(l) 

Endif 

write (6,*) ' Class threshold used is ',dgl 

!Read color image to be filtered: 

write (6, 150) inf il t ertv 

150 format (/' Input (dirty) file is ',a50) 

deltime-=Dtime (timearray) 

plml«malloc(iml bytes) 

If (piml .eq. OT then , , J , _ 

write (6,*)' Malloc for "iml" failed— abort' 
Call Exit(l) 

Endif 

Call Irgb load isl (infil, iml, nbx, ribands, npixm,nlinm, 
1 npix, nTin, pixformi, 0) 

pixformo«pixformi !Out datatype set equal to in type! 

write (6,160) clsfil 
160 formate Input class file is ',a50/> 
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pcls=malloc (els bytes) 
If (pels .eq. OT then 

write (6,*)' Malloc for "els" failed — abort' 

Call Exit(l) 

Endif 

Call Binput isl (clsfil, els, npixm, nlinm, npix, nlin, pixformc, 
l,nEc) 

pim2«malloc (im2 bytes) 
If (pirn? .eq. OT then 

write (6,*)' Malloc for "im2 w failed — abort' 

Call Exit (1) 

Endif 

pcls2=malloc (cls2_bytes) 
If (pcls2 .eq. 0) then 

write (6,*)' Malloc for "cls2" failed— abort' 

Call ExitU) 

Endif 

i 

inflg^l ! Start by assuming •els" is input class map! 

! Clean by erosion iteratively until no more pixels are replaced: 

write (6, *) ' [Beginning cleaning] ' 
cnttot=0 

iystart=l ! Required initializations 

iyend=nlin 

Do it«l,itmax 

cnt=0 

If (inflg .eq. 1) then 

Call Dedirt5x5 byte2_nopad (iml, els, im2, 

cls2, npixm, nlinm, npix, nlin, nbx, nbands, 
dgl , cnt , f min ) 
Call Dedirt5x5 byte2_nopad2 (iml, els, im2, 

cls7, npixm, nlinm, npix, nlin, nbx, nbands, 
2 dgl, cnt, fmin, iy start, iyend) 

inflg«2 

Else if (inflg .eq. 2) then 

Call Dedirt5x5_byte2_nopad(im2,cls2, iml, 

els, npixm, nlinm, npix, nlin, nbx, nbands, 
dgl, cnt, fmin) 
Call Dedirt5x5_byte2_nopad2 (im2, cls2, iml, 

1 els , npixm, nlinm, npix, nlin, nbx, nbands, 

2 dgl, cnt , fmin, iystart, iyend) 
inflg=l 

Else 

write (6,*)' Inflg out of range — abort' 
Call Exit(l) 

Endif 

cnttot=cnttot +cnt 
write (6,200) it, cnt 

formatC After it=',i3,' # pels cleaned = ',i6) 

If (cnt -eq. 0) then !No more dirt-class pels to be cleaned! 
go to 500 !Dreaded "go-to" statement! 

Endif 
Enddo 

write(6,*)' Warning-Cleaning loop exhausted after it=' , itmax 
write (6,*)' Aborting due to loop exhaustion' 
Call Exit(l) 
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500 continue 

write <6 f 550) it 
550 formate Cleaning ended after iteration^' , i3) 

write (6,552) cnttot 
552 formate Total # of points f iltered=' , 16) 

! Release classmap memory: 
Call Free (pels) 
Call Free(pcls2) 

If (inflg .eq. 1) then . 

Call Free(pim2) ! release unneeded image memory 
call Irgb_write_isl(outf il, iml, nbx,nbands, npixm r nlinm, 
npix, nlin, pixf ormo, 0) 
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Else 



Endif 



Call Free(piml) 'released unneeded image memory 

call Irgb_write_isl (outf ii, im2, nbx, nbands, npixm, nlinm, 

npix, nl in, pixf ormo, 0) 



deltiroe=Dtime (timearray) 
write (6, 800) timearray (1) , timearray (2) 
800 format-*/' Buster Erode_Am: : user time (sec):', 
1 f8.2, '; system time (sec) : ' , f 8 . 2/ ) 

write(6,8S0) outfil 
850 formatC Cleaned output image file is ',a50) 

end 
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Job: binput isl.vf 

Date: Wed Apr 21 15:08:40 1993 
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c Binput_isl: (10-7-92) Reads one band of an isl-format Image into 
c a single byte-type array. Input is clipped to range 0-255. 



70 



Subroutine Binput isl (inf il, im, npixm, nlinm, npix, nlin, pixf onn, 
band, nban<Js) 

Character*50 infil 
Character cval 

Byte im (npixm, nlinm) ,bval 

Integer*2 ibuf f (4096) , npix, nlin 

integer imgin, ostatus, gstatus 
integer opnirag,getdef , rdline 
integer nbands, pixf onn, band, iband, ival 
integer rstatus, npx, nln, npixm, nlinm 
integer arrtyp, lininc, jtype 



is Equivalence (cval,bval) 
i 



c include ' /local/include/iopackagef ' 

! Above comment-out must be removed after f77cvt conversion!!!! 

arrtyp=IDINT2 

20 c Open input image: 

imgin=0 

ostatus=opnimg (imgin, infil, 0, .false. ) 
If (ostatus .ne. SYSNRM) then 
write (6, 10) infil 
10 formate In image ',a50,' not open-abort') 

Call Exit(l) 

25 Endif 

gstatus^getdef (imgin, npx, nln, nbands, pixf orm) 
If (gstatus -ne. SYSNRM) then 
write (6, 20) 

20 formate Gstatus error on input — abort') 

Call Exit(l) 

30 Endif 

npix*npx 
nlin=nln 
write (6, 45) 

45 format {/' "Binput isl": Npix Nlin Type Nbands Band') 

write (6, 5V) npix, nlin # pixf orm, nbands, band 
50 format (14x,2i5,i6, iB, i5) 



35 



c.o. Read each line of input image: 

lininc«nlin/4 
jtype«=lininc 
rstatus«0 

40 Do j«=l,nlin 

c ... Read input lines & calc. sum: 

line»j-l 
iband=band -1 

rstatus-rdline (imgin, line, iband, ibuf f, 4096, arrtyp) 
If (rstatus .ne. SYSNRM) then 
45 write (7, 80) j 

80 formate Error reading line -abort') 

Call Exit (1) 
Endif 



Do i*l # npix 
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Enddo 



ival-MaxO (ibuff (i) , 0) 
ival^MinO (ival, 255) 
cval=Char (ival) 
ira(i, j)=bval 



Enddo 



!End line loop! 



c Close images 

Call clsimg (imgin) 

return 

end 
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Job: dedirt5x5 byte2 nopad dugas.vf 
Date: Wed Apr 21 15:07:51 1"?93 
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c Dedirt5x5_byte2_nopaa: 9 

c Eliminates "dead zone** 

c of 2 pixels around image edges (i.e. cleaning now takes place at all 

c pixels of the image) . Only need one clean pixel in a 5x5 neighborhood t 

c clean a dirty border pixel (i.e. nrain edge^l) . 

5 

c Output non-dirt or cleaned classes in classmap have a value of "0". 

c Note that within one call of Dedirt, cleaned GLs are *not* used 

c to clean their neighbors; this will, however, occur in 

c subsequent iterations (further calls to dedirt) . 

c Ad hoc "dirt" remover. Assumes that dirt is detectable as 

70 c a simple grey-level threshold in the input classmap. 

c Filters all pixels *above* a certain 

c threshold grey level ("gth") . 

c Note that a 1-pixel dilation of the classmap is suggested prior 

c to first using this subroutine, to avoid propagation of dirt 

c "ghosts'*. 

, 5 C 

c Image arrays: (Max dims :npixm, nlinm; Actual dims: npix,nlin) 

c iral : Input (dirty) 3-band image 

c els: Input dirt 1-band classmap 

c cls2: Output dirt 1-band classmap 

c im2: Output (cleaned) 3-band image 

c 

c Input: np,nl :Size of filter window. 

20 c dgl : Threshold classmap GL value req' d for cleaning. 

c fmin rFraction of window needed to be "clean" to be filtered 

c (l=<nmin<np* (nl-1) /2 ). 

c Output :cnt :No. of pixels which were filtered. 

Subroutine Dedirt 5x5 Byte2_Nopad (iral, els, im2, cls2, npixm, nlinm, npix, nlir. 
25 1 dgl,cnt, fminT 

Integer*2 iml (nbx, npixm, nlinm) , im2 (nbx, npixm, nlinm) 

Integer*2 npix,nlin 

Integer npixm, nlinm, np,nl,np2,nl 2, cnt,dcnt, iO, ix, jO, jx,dgl 
Integer nbx, nb, num_neigh, nmin, delx, dely 
30 Integer nmin_edge,runin_cent 

Byte els (npixm, nlinm) , cls2 (npixm, nlinm) ,bdgl 

Real*4 sum (3) ,wtsum, wt (3, 3) , fmin 

! fmin= .24 implies nmin=6 for 5x5 window. 

35 Character cdgl 

Data np, nl/5,5/ ! Erosion window size! 

Data vt(l,l),vt(2,l),wt (3, 1) /l . , 1. , 0.5/ !Wts specific t: 

Data wt (1,2), wt(2, 2), wt (3, 2)/l., 0.707, 0.4472/ 
Data wt (1,3) ,wt (2,3) ,wt (3,3) /0. 5, 0.4472, 0.3536/ 
Data nmin_edge/0/ 
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Equivalence (cdgl,bdgl) 



cnt=0 
dcnt-0 
nl2«nl/2 
np2~np/2 

nmin_cent=Nint (fmin*np*nl) 
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cdgl=Char (dgl) 10 
write (6, 10) 

formate Dedirt5x5 Byte2_Nopad: Npix Nlin Np Nl Fmin Nmin Dgl') 
write {6, 11) npix, nTin, np, nl, fmin, nmin cent, dgl 
f onnat <23x,2i5, 2i4, f 8 . 3, i5, i6/) 



! Slide window across image: 
Do j-l,nlin 
Do 1-1, npix 

!i, j are coordinates of center of window! 

•Test local value of classmap: 

If (cls(i,j) .eq. bdgl) then 'Filter the point 
dcnt=dcnt +1 
j0=j-nl2 
jx=j+n!2 
jO=MaxO(jO,l) 
x=MinO < jx,nlin) 
0~i-np2 
ix«i+np2 
iO=MaxO<10, 1) 
ix«MinO (ix, npix) 

!iO,ix, j0,jx are coordinates of extrema of 
! current window 
Do ib=l,nb 

sum(ib) =0 . 

Enddo 
wtsum=0- 
num_neigh=0 

!Calc. over local ' , dirty ,, window: 
Do j2=j0,jx 

Do i2=i0,ix 

!If neighbor is *not* dirt, then: 
If <cls(i2,j2) -ne. bdgl) then 

!Calc. a distance-weight for that GL: 
!dis= (Float (i -i2>)**2 + (Float (j -j2))**2 
!dis=Sqrt (dis) 
!wt-l./dis 
! wtsum=vrtsum +wt 

delx=Int (Abs(i -i2) +1.) 
dely=lnt (Abs(j -j2) +1.) 
wtsum=wtsum +wt (delx, dely ) 

Do ib«l,nb 

sum (ib) »sum (ib) +wt (delx, dely) *ial (ib, 

Enddo 

! Count no. of non-dirt window pixels: 
num neigh=num neigh +1 
Endif 

Enddo 
Enddo 

!Calc. local value of ■nmin": 

If (j .le. nl2 .or. j .ge. (nlin -nl2 +1)) then 

nmin ■nmin edge !Top or bottom border! 
Else if (i ,le. np2 -or. i -ge. (npix -np2 +1)) then 

nmin«nmin- edge !Left or right border! 
Else " 

nmin*nmin_cent !Not in a border region! 
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adif 11 

!If have mln. no. of non-dirt neighbors, then 
! replace dirty GL w/ wted ave of nondirty GLs: 
If (nun neigh .gt. nmin) then 'Replace the center GL: 
5 Do i5=l , nb 

im2 <ib, i, j) =Nint (sure <ib) /wtsum) 

Enddo 

cls2(i,j)=0 !Don't filter point again! 
cnt»cnt +1 

Else !Not enough non-dirt pixels within window! 
io !Keep the old (dirty) GL for now: 

Do ib=l, nb 

im2 <ib. i, j)=iml (ib, i, j) 

Enddo 

cls2(i, j)=cls(i, j) 
Endif 

?5 Else 'Point is non-dirt , no interpolation! 

•Copy input GLs to output: 
Do ib«l,nb 

im2Ub,i, j)=iml(ib,i, j) 

Enddo 

cls2(i,j>«cls(i, j) 

- Endif 

20 Enddo !End "i" loop 

Enddo !End "j" loop 

write (6,500) dcnt,cnt 
500 formate No. of dirt pels:',i6,'; no. of cleaned pels:',i6) 



25 



If (dent .gt. 0 -and. cnt .eq. 0) then 

write (€,*>' Warning — stable w/ some dirt remaxning ! ! ! 
write (6,*)' [aborting]' 
Call Exit(l) 

Endif 

return 

end 
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c li*gb_vrite_isl : Writes a 3-band image from a Int*2 array. 

Subroutine Irgb write is 1 (bname,b, nbx, nb, npixm, nlinm, npix, nlin, 
1 pixf orm, lunprl:) 

Character*50 bname, outf il 

Integer* 2 b (nbx, npixm, nlinm) , npix, nlin, ibuf f (4 096, 3) 

integer imgout, pstatus, ©status 

integer ribands, pixf orm, nbx, nb, npixm, nlinm 

integer wstatus, npx, nln, ival 

integer arrtyp, lininc, jtype,pixdira, lunprt 

integer putdef , opnimg, wrline 

Data pixdim/4096/ 

include ' /local/include/ iopacxagef 

arrtyp=IDINT2 
nbands=nb 

write(6,45) npix, nlin, nbands, pixf orm 
45 format <2x, ' w Irgb_write_isl- : : Npix: ',15,'; NIin:',i5,'; Nb:',i2, 

1 '; Type:',i2) 

If (lunprt .gt. 0) then 

write (lunprt, 45) npix, nlin, ribands, pixf orm 

Endif 



c . . - Open output image: 

imgout=0 
npx=npix 
nln=nlin 
outf il=bname 

25 pstatus=putdef (imgout , npx, nln, nbands , pixf orm) 

If (pstatus .ne. SYSNRM) then 

write (6,*)' Pstatus error for output — abort' 
stop 

Endif 

ostatus=opnimg (imgout, outf il, 1, -false. ) !"1" for write-only! 
If (ostatus .ne. SYSNRM) then 
30 write(6,*)' Output image not open — abort' 

stop 

Endif 

c... Write each line of output image: 

lininc=nlin/5 
jtype=lininc 
35 Do j=l,nlin 

c ... load 2-d array & write lines: 

Do i=l,npix 

Do ib=l, nbands 

ibuf f (i, ib) =b (ib, i, j) 

Enddo 
Enddo 

If (pixform .eq. IDBYTE) then 
Do ib*=l, nbands 
Do 1=1, npix 

ival=MaxO (ibuff (i,ib) ,0) 
ibuf f ( i , ib ) ~MinO ( ival , 255) 

Enddo 
Enddo 

Else if (pixform .eq. IDINT1) then 
Do ib*=l, nbands 
Do i=l,npix 

ival-MaxO (ibuf f (i, ib) , -128) 
ibuff (i,ib)-MinO(ival, 127) 

Enddo 

50 Enddo 
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Else if (pixrorm .eq. IDINT2) th^n 
continue 

Else if (pixform .eq. IDREA4) then 
continue 

Else 

write (6, 100) pixform 
100 formate Pixformo=' ,i3, ' not implemented — stop ) 

stop 

Endif 

line=j -1 

wstatus=wrline (imgout , line, -1, ibuf f , pixdim, arrtyp) 
If (wstatus .ne. SYSNRM) then 

write (6,*)' Error in wstatus — abort' 

stop 

Endif 

c - If < j -ge. jtype) then 

c write(6,200) j 

c 200 format <10x, ' Line has been written') 

c jtype=jtype +lininc 

c Endif 

Enddo 'End line loop! 

c. Close image: 

Call clsimg (imgout) 

return 

end 
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Claims 

1. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
30 pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

as d) forming a residual value for each pixel as a function of the difference between the corresponding 

feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
. residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
40 b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

2. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
45 pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

so d) forming a residuals image which consists of the difference between the feature image pixel values 

and the filtered values of step c); 

e) testing each pixel value of the residuals image of step d) to determine if the residual value is 
within a range of residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
55 b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

3. A method for the detection and removal of defects in digital images comprising the steps of: 
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a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

The method according to Claim 1, 2, or 3. wherein the map formed in step f) is conditioned by 
morphological dilation of the pixels marked as defect. 

The method according to Claim 1, 2, or 3, wherein the map formed in step f) is conditioned by 
morphological erosion of the pixels not marked as defect. 

The method according to Claim 1 , 2, or 3, wherein the map of step f) is formed with each pixel marked 
according to the magnitude of its values according to steps b) and e). 

A method for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising the steps of: 

a) preprocessing a lower spatial resolution representation of the digital image by forming a feature 
image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; 

g) correcting the lower resolution representation of the digital image using the map created in step 
f);and 

h) correcting the higher resolution representations of the digital image using the map created in step 
0. 

An apparatus for the detection and removal of defects in digital images comprising: 

means for preprocessing a digital image represented by pixel values to form a feature image in 
which defect pixels have large local contrast from their neighboring non-defect pixels; 

first means for testing the value of each feature image pixel value to determine if the value is within 
a range of pixel values expected of defects; 

filter means for edge-preserving spatial filtering of each of the feature image pixel values; 

means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the filtered value from said filter means; 

second means for testing each residual value from said forming means to determine if the residual 
value is within a range of residuals values expected of defects; 

means for forming a map in which each pixel that is within the range of expected defects according 
to said first and said second testing means is marked as a defect; and 

means for correcting the digital image as a function of the formed map. 

An apparatus for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising: 

means for preprocessing a lower spatial resolution representation of the digital image by forming a 
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feature image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

first testing means for testing the value of each feature image pixel value to determine if the value 
is within a range of pixel values expected of defects; 

filter means for edge-preserving spatial filtering of each of the feature image pixel values that are 
within the range of expected defects; 

residual means for forming a residual value for each pixel as a function of the difference between 
the corresponding feature image pixel value and the filtered value of said filter means; 

second testing means for testing each residual value from said residual means to determine if the 
residual value is within a range of residuals values expected of defects; 

mapping means for forming a map in which each pixel that is within the range of expected defects 
from said first and said second testing means is marked as a defect; 

means for correcting the lower resolution representation of the digital image using the map created 
by said mapping means; and 

means for correcting the higher resolution representations of the digital image using the map 
created by said mapping means. 

10. A method for the detection and removal of local defects in digital images comprising the steps of: 

a) generating a feature image from an original digital image; 

b) creating an EPS-residuals image using the feature image; 

c) testing and mapping image pixels based upon the values of the EPS residuals image and the 
pixel values of the feature image; 

d) modifying the mapped image as a function of mapped values of neighboring pixels; 

e) cleaning of defect pixels by an EPS-substitution process or by a greyscale erosion process; and 

f) cleaning of defect pixels in different spatial resolution representations of the original image if 
different spatial resolution representations are present. 

11. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image to form a feature image; and 

b) testing each feature image to determine if a value is within a range of pixel values expected of 
defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects. ' 
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© A technique for the detection and removal of local defects in digital continuous-tone images. 

(omen 



® The present invention is a method for automati- 
cally detecting and correcting a wide range of local 
digital image defects with minimal user intervention. 
The detection process employs brightness and color 
thresholds in conjunction with magnitude thresholds 
on residuals of nonlinear spatial filters to separate 
defects from scene content with minimal confusion. 
The detected defects are then cosmetically correct- 
ed by combinations of nonlinear smoothing and 
grey-scale erosion. Several options are outlined for 
the feature selection, detection, and cleaning oper- 
ations depending on source type and computational 
constraints. 
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